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ABSTRACT  An  efficient  finite  element  procedure  to  analyze  transient  phenomena 
in  dry  and/or  fluid-saturated  porous  media  is  presented.  The  saturated  porous  medium 
is  modeled  as  a  twO-phase  system  consisting  of  a  solid  and  a  fluid  phase.  Time  integration 
of  the  resulting  semidiscrete  finite  element  equations  is  performed  by  an  implicit-explicit 
algorithm.  In  order  to  remove  the  time  step  size  restriction  associated  with  the  presence  of 
the  stiff  fluid  in  the  mixture,  the  fluid  contribution  to  the  equations  of  motion  is  always 
treated  implicitly.  The  procedure  allows  an  optimal  selection  of  the  time  step  size 
independently  of  the  fluid.  Depending  upon  the  particular  intended  applications  (e.g., 
seismic,  blast  loading),  the  fluid  may  be  assumed  incompressible  or  compressible. 

Accuracy  and  versatility  of  the  proposed  procedure  are  demonstrated  by  applying  it 
to  analyze  a  number  of  dynamic  soil  and/or  soil  structure  interaction  tests  and  quasi-static 
consolidation  tests  performed  in  centrifuges. 
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INTRODUCTION  AND  PROBLEM 


The  Navy  has  $25  billion  worth  of  facilities  in  seismically  active 
regions.  Each  year  $200  million  of  new  facilities  are  added  to  those  in 
seismically  active  areas.  The  Navy,  because  of  its  mission,  must  locate 
at  the  waterfront  with  a  high  watertable  and  often  on  marginal  land. 
Seismically  induced  liquefaction  is  a  major  threat  to  the  Navy.  Pres¬ 
ently,  procedures  do  not  exist  to  analyze  the  effect  of  liquefaction  on 
structures.  Developing  an  effective  stress  soil  model  will  provide  a 
tool  for  such  analysis  of  waterfront  structures. 

To  understand  the  significance  of  liquefaction,  it  is  important  to 
note  the  damage  caused  in  recent  earthquakes.  The  following  summarizes 
recent  experiences  during  major  earthquakes. 

1960  Chilean  Earthquake  (Magnitude  6-8.3) 

Most  spectacular  damage  occurred  in  Puerto  Montt,  to  quay  walls, 
steel  sheet  piles,  and  sea  walls.  Liquefaction  of  the  loose  fine 
sandy  soils  was  the  primary  cause  of  the  failures. 

1964  Alaska  Earthquake  (Magnitude  8.4) 

Severe  damage  at  Anchorage,  Cordova,  Homer,  Kodiak,  Seldovia, 
Seward,  Valdez,  Klawock,  and  Whittier.  Large-scale  land  slides  and 
liquefaction  induced  most  of  the  extremely  heavy  damage  and  total 
destruction. 
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Severe  damage  in  Niigata  Port  (West  Harbor).  Areas  affected  were 
Additional  Harbor,  Yamanoshita  Wharf,  North  Wharf,  East  Wharf, 
Central  Wharf,  South  Wharf,  Kurinoki  River  Landings,  Bandai  Island 
Wharf,  Shinano  River  Left  Bank  Bulkhead,  and  West  Coast  Bulkheads. 


Liquefaction  caused  most  of  the  heavy  damage- 


1968  Tokachi-Qki  Earthquake  (Magnitude  7.8) 


Ports  affected  were  Hachinohe,  Aomori,  Hakodate,  and  Muroran. 
Damage  was  relatively  light  compared  to  that  caused  by  Niigata 
Earthquake.  Most  of  the  damage  occurred  to  structures  of  smaller 
scale.  Liquefaction  was  not  the  primary  cause  of  damage  even 
though  spouting  sand  sediments  were  seen  at  several  waterfront 
areas  near  the  damaged  structures. 


1973  Nemuro-Hanto-Qki  Earthquake  (Magnitude  7.4) 


Severe  damage  occurred  mainly  in  Hanasaki  and  Kiritappu  Ports. 
Nemuro  Port  situated  only  6  km  away  from  Hanasaki  Port  sustained 
very  slight  damage.  The  damage  was  attributed  to  soil  lique- 
faction. 


1978  Mi vagi-Ken-Qki  Earthquake  (Magnitude  7.4) 


Areas  affected  were  Shiogama,  Sendai,  and  Ishinomaki  Ports,  and 
Ishinomaki  and  Yuriage  Fishing  Ports.  The  damage  in  Ishinomaki 
Port  accounted  for  approximately  90  percent  of  the  total  damage 
costs  at  port  and  harbor  facilities  caused  by  this  earthquake. 
Gravity  quay  walls  and  piers  suffered  various  degrees  of  damage. 
Sheet  pile  quay  walls  damaged  were  primarily  due  to  liquefaction 
of  fill  materials.  Liquefaction  again  played  a  significant  role 
in  this  earthquake.  At  sites  where  liquefaction  did  occur,  the 


damage  to  the  port  and  harbor  structures  was  very  severe. 
Conversely,  the  damage  to  port  and  harbor  structures  was  small  at 
sites  where  no  liquefaction  occurred. 

As  can  be  seen,  liquefaction  played  a  major  role  in  waterfront 
damage,  most  of  the  time  being  the  single  cause  of  widespread  losses. 
Fortunately  the  United  States  has  not  suffered  a  devastating  earthquake 
in  recent  years.  However,  the  seismic  risk  is  great,  particularly  in 
the  West  where  it  is  estimated  that  there  is  a  5  percent  annual  prob¬ 
ability  of  a  major  event  in  Southern  California  that  could  affect  a 
number  of  Naval  bases.  This  problem  has  been  noted  in  an  ONR  sponsored 
study  evaluating  the  Navy's  seismic  vulnerability.  The  experience  noted 
in  recent  earthquakes  is  that  liquefaction  greatly  increased  the  amount 
of  damage  observed  in  waterfront  facilities.  Particular  problems  exist 
with  sheet  piles,  quay  walls,  wharfs,  and  embedded  structures.  Conven¬ 
tional  buildings  also  suffer  severe  damage. 

The  Navy's  experience  has  been  limited  to  damage  inflicted  in  the 
1964  Alaskan  earthquake.,  heavy  damage  was  noted  in  the  seawall  at  the 
Kodiak  Naval  Station.  One  foot  of  differential  settlement  was  noted 
beneath  aircraft  hangars.  Compaction  of  fill  occurred  under  asphalt 
aircraft  ramps.  It  is  significant  to  note  that  these  facilities  were 
constructed  on  15  to  20  feet  of  engineered  fill  where  seismically 
induced  pore  pressure  increases  would  be  expected  to  reduce  soil  stiff¬ 
ness  and  shear  strength.  The  damage  noted  was  caused  by  soil  failure, 
and  in  addition,  substantial  damage  was  caused  by  the  seismic  sea  wave. 
The  United  States  has  not  had  a  large  number  of  events  exposing  Navy 
facilities  to  damage.  However,  the  Japanese  have  had  a  number  of  events 
and  their  experience  illustrates  that  seismic  liquefaction  was  respon¬ 
sible  for  most  waterfront  damage. 

The  conclusions  from  this  are: 

(1)  Seismic  liquefaction  causes  severe  damage  to  waterfront 
structures . 
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(2)  The  Navy,  located  in  seismically  active  areas  having  numerous 
waterfront  structures  on  marginal  land,  is  vulnerable  to  substantial 
damage. 

(3)  Techniques  presently  do  not  exist  to  accurately  analyze  the 
response  of  a  large  complex  waterfront  structure  on  soil  in  which  seis¬ 
mically  induced  pore  pressures  cause  loss  of  soil  stiffness  and  shear 
strength  (liquefaction). 

The  most  promising  solution  to  this  problem  is  developing  a 
constitutive  soil  relationship  that  is  capable  of  accurately  predicting 
soil  behavior  under  generalized  loading  conditions.  Implementing  this 
effective  stress  soil  model  into  a  finite  element  computer  program  would 
allow  analysis  of  soil  and  structure  together. 

TECHNICAL  BACKGROUND 

The  analysis  of  dynamic  transient  phenomena  in  fluid-saturated 
porous  media  is  of  great  interest  in  geophysics  and  geotechnical 
engineering.  Fluid  saturation  of  an  otherwise  inviscid  porous  solid 
skeleton  introduces  a  time  dependence  into  the  reponse  to  applied  loads. 
Biot  (Ref  1)  first  considered  the  propagation  of  harmonic  waves  in  a 
fluid-saturated  porous  medium.  Since  then,  his  theory  and  results  have 
been  the  standard  reference  and  basis  for  most  of  the  investigations  in 
acoustics,  geophysics,  and  geomechanics.  Many  one-dimensional  wave 
propagation  theories  have  since  been  proposed  (see  Ref  2  and  3)  for 
recent  surveys  of  western  and  Russian  literature),  and  one-dimensional 
wave  propagation  numerical  results  were  first  presented  in  References  4 
and  5.  The  need  for  a  general  multidimensional  formulation  and  solution 
technique  has  become  important  in  recent  years  because  of  the  increased 
concern  with  the  dynamic  behavior  of  saturated  soil  deposits  and  associ¬ 
ated  liquefaction  of  saturated  sand  deposits  (see  Ref  6  and  7)  under 
seismic  loading  conditions.  Also  concern  in  marine  foundation  engineer¬ 
ing  with  water  wave  induced  dynamic  pore  pressures  in  saturated  marine 
deposits  has  spurred  interest  in  the  subject  matter  (see  Ref  8  and  9 


for  related  analytical  solutions).  Most  of  the  solution  procedures 
reported  in  the  literature  are  restricted  to  linear  systems  Ghaboussi 
and  Wilson  (Ref  10)  first  proposed  a  multidimensional  finite  element 
numerical  scheme  to  solve  the  linear  coupled  governing  equators 
Despite  the  extensive  literature  published  in  soil  dynamics  (see 
Reference  11  for  extensive  references)  no  general  technique  capable 
of  accounting  for  all  present  nonlinear  effects  (large  deformations' 
strains,  nonlinear  material  behavior)  has  yet  been  fully  developed  and 
implemented,  although  attempts  at  presenting  a  suitable  general  frame 
work  have  been  reported  (see  Ref  12,  13,  14,  16,  16,  and  17) 

Professor  J.H.  Prevost  of  Princeton  Jniversity  (Ref  13.)  has 
developed  an  efficient  finite  element  numerical  procedure  to  analyze 
transient  phenomena  in  fluid  saturated  porous  media  The  saturated 
porous  medium  is  modeled  as  a  two-phase  system  consisting  of  a  solid  and 
a  fluid  phase.  The  solid  skeleton  may  be  linear,  or  nonlinear  and 
hysteretic.  Large  deformations  may  also  be  included  The  fluid  may  be 
compressible  or  incompressible  depending  upon  the  intended  applications 
(eg.,  seismic,  blast  loading).  Time  integration  of  the  resulting  semi¬ 
discrete  finite  element  equations  is  performed  by  using  an  implicit- 
explicit  algorithm  (Ref  18  and  19).  In  order  to  remove  the  time  step 
size  restriction  associated  with  the  presence  of  the  stiff  fluid  in  the 
mixture,  the  fluid  contribution  is  always  treated  implicitly. 

This  study  is  directed  toward  examining  the  predictive  capabilities 
of  the  numerical  procedure  proposed  in  Reference  20.  °f  particular 
interest  is  the  validity  of  the  proposed  numerical  model  in  capturing 
adequately  the  generation  and  dissipation  of  excess  pore-water  pressures 
in  saturated  sand  deposits  during  (and  after)  earthquakes,  and  its  per¬ 
formance  in  dynamic  soi 1  -  structure  interaction  problems.  The  most 
appropriate  method  for  such  a  validation  study  would  be  to  utilize  field 
data  from  instrumented  prototype  situations.  However  such  a  study  is 
preempted  by  the  paucity  and  scarcity  of  the  field  data.  In  the  absence 
of  actual  prototype  earthquake  field  data,  an  alternate  method  of  vali¬ 
dation  is  provided  by  analyzing  centrifuge  soil  model  test  data. 
Although  imperfect  in  many  respects,  it  is  felt  that  dynamic  centrifuge 
soil  model  tests  can  still  provide  a  data  base  for  calibration  of 
numerical  procedures.  A  number  of  dynamic  centrifuge  soil  model  tests 
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interaction  of  two  such  phases  The  required  formalism  for  the  opvp  of.  ■ 

ment  of  general  nonlinear  equations  was  later  introduced  through  the 
theory  of  mixtures  (Ref  26,  27,  and  28)  General  mixture  results  can  be 
shown  through  formal  1 i nea r i za t i on  of  the  field  and  constitutive  equa¬ 
tions,  to  reduce  to  Biot  linear  poroelastic  model  (see  Ref  29).  The 

balance  laws  for  the  two-phase  mixture  are  summarized  in  this  report 


Balance  of  Mass 


On  or  _  ..a 

-zr-  -*■  n  V  •  V 


(n  )  (a  =  s.w) 


Balance  of  Linear  Momentum 


„  a  ^  M  a  .  a  a  ,  .  ... 

V-a+g+pb  =  pa  (a  =  s,w)  (2) 

a 

where:  D/Dt  =  material  derivative  following  the  motion  of  the 
a-phase 

ap  =  napa  =  macroscopic  average  mass  density 
Pa  =  microscopic  mass  density 

na  =  fraction  of  elemental  volume,  dV,  occupied  by  a-phase 
(i.e.  ,  na  =  dV°7dV). 

Clearly,  I  n°  =  ns  +  nW  =  1, 

where:  s  and  w  =  solid  and  fluid  phases,  respectively 
nw  =  porosity 

av  =  velocity  (spatial)  of  a-phase 
aa  =  acceleration  of  a-phase 
b  =  body  force  per  unit  mass 

^  Ql 

g  =  momentum  supply  to  the  a-phase  from  the  other  phase, 
subject  to  I  g  =  gs  +  gw  =  0. 

In  the  following,  momentum  interaction  consists  of  diffusive  and  dilata- 
tional  contributions,  viz. , 

E  =  -  4  .  (v  -  v  )  -  pw  V  n  (3) 

where:  %  -  symmetric,  positive  definite  second-order  tensor 

Pw  =  fluid  pressure 

The  first  term  is  sometimes  called  the  "Stokes  drag"  (see  Ref  26). 

In  Equation  2,  aa  =  partial  (Cauchy)  stress  tensor  corresponding  to 

the  a-phase.  The  partial  stress  tensor  oW  corresponding  to  the  fluid 

w 

phase  is  equal  to  n  times  the  pore  fluid  stress  ctw,  i.e.. 


(4) 


w 


w 

=  n  a 
~w 


However,  the  partial  stress  tensor  a  corresponding  to  the  solid  phase 
is  not  the  effective  stress  o's  of  classical  soil  mechanics  (Ref  25), 
but  rather  is 


i  s  s 
=  o  +  n  o 

~  ~w 


(5) 


o  accounts  for  the  effects  of  the  pore  fluid  stress  on  the 
~w 


where  n 

individual  solid  grains  which  constitute  the  solid  skeleton.  The  global 
stress,  a,  is  the  sum  of  the  partial  stresses,  and  is  equal  to  (from 
Equations  4  and  5) 


o  =  o  +o+o  +o 


(6) 


as  postulated  in  classical  soil  mechanics  (Ref  25). 


Equations  of  State 


Constitutive  Assumptions 


A  rate-type  constitutive  equation  describes  the  behavior  of  the 
porous  solid  skeleton,  of  the  following  form: 


s 

_D_ 

Dt 


(a'S)  =  DS  :  v*  ^  +  DG  :  v^  ^  -  (o,S  x  6's)  :  ^ 


(9a) 


where:  v 


(  )  =  symmetric  parts  of  the  solid  velocity  gradient 
Vj-Sj  =  skew- symmetric  parts  of  the  solid  velocity  gradient 
Ds 

~  =  material  constitutive  tensor,  an  objective  tensor 
valued  function  of  possible  o'  and  the  solid 
deformation  gradient 


Namely, 


=  contribution  from  the  rotational  component  of  the 
stress  rate. 


Di  jk£  2Cct i£  6jk  +  aj£  6ik  '  aik  6j £  "  ajk  6i£^ 


..s 


•  s 


(9b) 


The  last  term  on  the  right-hand  side  in  Equation  9a  is  introduced  to 
ensure  that  the  tangent  stiffness  operator  obtained  through  lineariza¬ 
tion  of  the  momentum  equations  possesses  the  major  symmetry  as  the  Ds 
tensor.  Many  nonlinear  material  models  of  interest  can  be  put  in  the 
above  form  (e.g. ,  all  nonlinear  elastic  and  many  elastic-plastic 
material  models).  The  particular  form  of  the  constitutive  equation 
(Equation  13)  adopted  here  was  first  proposed  by  Hill  (Ref  30)  within 
the  context  of  plasticity  theory.  Appropriate  expressions  for  the 
effective  modulus  tensor,  0s,  for  soil  media  are  given  in  References  31 
and  32.  For  a  linear  isotropic  elastic  porous  skeleton: 


Dijk£ 


Xs  6.  . 


U  6k£  +  V  (6ik  6j£  +  6i£  V 


(10) 
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where:  \S ,  VS  =  effective  Lame's  moduli 


6.  .  =  Kronecker  delta 
ij 

The  following  constitutive  equation  is  assumed  to  describe  the 
behavior  of  the  fluid  phase 

a  =  -  p  6  (11) 

where:  p^  =  pore-fluid  pressure 

We  assumed  that  the  fluid  has  no  average  shear  viscosity.  Further,  the 

fluid  flow  is  assumed  barotropic  so  that  the  fluid  kinetic  equation  of 

state  is  independent  of  the  temperature,  viz. , 

F(PW.PW)  =  0  (12) 

from  which  it  follows  that 


tw(0w>  =  ^5  (13) 

W  A 

where:  Aw  =  PW^PW/3PW  =  bulk  modulus  of  the  fluid  phase. 

The  fluid  pressure  can  thus  be  determined  from  Equation  8,  which  now 
writes: 

w 

Dt  (pw)  =  -^7  &  •  vW)  +  V  •  C(1  -  nW)  Vs]}  (14) 

1  v 

For  soil  media,  the  compressibility  of  the  fluid  phase  is  often  no 
smaller  than  the  compressibility  of  the  solid  skeleton.  Therefore,  the 
fluid  phase  may,  in  some  soil  dynamics  applications,  be  regarded  as 
incompressible,  and  Equation  8  reduces  to: 


V  •  [nw  vw]  +  V  •  [(1  -  nw)  vs]  =  0 


Under  the  assumptions  described  above,  the  linear  momentum  equa¬ 
tions  (Equation  2)  simplify  to. 

PS  £S  =  £  •  2|S  -  nS  V  pw  -  ^  •  (vS  -  vW)  +  pS  b  (16a) 

s 

p  Dt  )  p(-  "  l  )  *  1  1  'nHPw  (16b) 

+  C  •  (v  -  V  )  +  p  b 

when  the  movement  of  the  solid  phase  is  used  as  the  reference  motion. 
When  inertia  and  convective  terms  are  neglected,  Equation  16b  reduces  to 
Darcy' s  law  as: 

nW(vW  ~  vS)  =  ~  (nW)2  il  •  (V  p  -  p  b)  (17) 

w  2  1 

and  thus  k  =  (n  )  yw  £  =  Darcy  permeability  tensor  (symmetric,  positive 

definite),  (units.  L/T),  yw  =  gp^  =  unit  weight  of  the  fluid; 
g  =  acceleration  of  gravity. 

Weak  Form  -  Semi -Discrete  Finite  Element  Equations 

The  initial  boundary  value  problem  consists  of  finding  the  solid 
s  w 

displacement,  v  ,  the  fluid  velocity,  v>  ,  and  the  fluid  pressure,  pw, 
(all  functions  of  position  and  time)  satisfying  the  field  equations 
(Equation  16a  and  16b)  together  with  the  constitutive  relations  and 
continuity  conditions  subject  to  appropriate  initial  and  boundary 
conditions.  In  order  to  reduce  the  number  of  unknowns,  the  fluid 
pressure  is  eliminated  from  the  formulation,  thereby  producing  a  most 
efficient  scheme.  In  the  case  of  a  compressible  fluid,  the  fluid 
pressure  is  determined  from  the  computed  velocities  through  time  inte¬ 
gration  of  Equation  14.  In  the  case  of  an  incompressible  fluid,  a 
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penalty-function  formulation  of  the  continuity  constraint  expressed  by 
the  storage  equation  (Equation  15)  is  used  to  compute  the  fluid  pressure 
as: 


w 

Pw  =  ■  ^  fX  *  +  l  •  [(1  -  Vs]}  (18) 

n 

w 

where  A  >  0  is  a  penalty  parameter,  not  the  effective  bulk  modulus 
appearing  in  Equation  14.  The  penalty  parameter  is  selected  as  a  large 
number.  This  parameter  is  further  discussed  later. 

The  weak  formulation  associated  with  the  initial  boundary  value 
problem  is  obtained  by  proceeding  along  standard  lines  (see  Ref  33). 
The  associated  matrix  problem  is  obtained  by  discretizing  the  domain 
occupied  by  the  porous  medium  into  nonoverlapping  finite  elements. 
Associated  with  this  discretization  are  nodal  points  at  which  shape 
functions  are  prescribed.  Two  sets  of  shape  functions  will  be  required 
for  the  solid  displacement  and  the  fluid  velocity  fields,  respectively. 
However,  since  attention  in  the  following  is  restricted  to  low  order 
(i.e.,  four-node  plane.,  eight-node  brick)  finite  elements,  which  are 
the  most  efficient  in  nonlinear  analysis,  the  same  shape  functions  are 
used  for  both  the  solid  and  the  fluid.  The  shape  functions  for  the 

solid  displacement  and  fluid  velocity  associated  with  node  A  are 

A  A  B  B 

denoted  by  N  .  They  satisfy  the  relations  N  (x  )  =  in  which  x 

d  notes  the  position  vector  of  node  B,  6^  =  Kronecker  delta.  The 
solution  of  the  Galerkin  counterpart  of  the  weak  formulation  is  then 
expressed  in  terms  of  the  shape  functions  and  gives  rise  to  the  follow¬ 
ing  system  of  equations 


MS  0 
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z  -  z 
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FW 
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V 
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where:  M0  =  mass  matrix 

a°  =  acceleration  vector 


v  =  velocity  vector 
Fa  =  force  vector 

Several  computational  simplifications  result  in  using  a  diagonal  mass 
matrix,  and  a  "lumped"  mass  matrix  is  used  throughout.  For  the  two- 
dimensional  (three-dimensional),  four-  (eight-)  node  bilinear  (tri- 
1  inear)  isoparametric  element, 

(M^j)  =  6AB  J  p“  NA  dO  (no  sum  on  A)  (20) 


where:  M .  .  =  the  elemental  mass  contribution  to  node  A  from  node  B  for 

1J  e 

directions  i  and  j  to  the  global  mass  matrix  0  =  spatial  domain  occu¬ 

pied  by  element  e. 

In  Equation  19,  Z  is  a  damping  matrix  arising  from  the  momentum 
transfer  terms  in  Equations  16  and  17  as: 


ZAB  =  f  NA.£  NB  do 

’J  J  ir 


The  solid  force  vector  F  is: 


(  Fext)  -  NS 


gxt  s 

where:  (F  )  =  vector  of  external  solid  forces  (i.e.,  body  force, 

surface  tractions) 


I 
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N  =  vector  of  solid  stress  forces,  viz., 


NA  5 


=  /  NA  .(o'. *  -  nS  p  6.  .) 

J  .  J  1 J  "  ’J 


The  fluid  force  vector  F  is: 


'  FeXt  J  -  NW 


-  ■  T*V'  ■_  ■, 


*  A  .%  ■«  a 


m 


Vet' 


0X  t.  w 

where:  (F  )  =  vector  of  external  fluid  forces 


N  =  vector  of  convective  and  fluid  stress  forces,  viz. , 


NA  )  =  f  pw  NA(v"  -  Vs!)  v?  .  dfi  -  [nA.  nW  p  dfl 

W  J  J  J  i  ,J  J  ,1  w 


Time  Integration 


Time  integration  of  the  semidiscrete  finite  element  equations 
(Equation  19)  is  performed  by  using  a  finite  difference  time  stepping 
algorithum.  Many  types  of  time  stepping  algorithms  and  algorithmic 
strategies  are  presently  available  (see  Ref  34  for  a  description  of 
most  widely  used  computational  transient  analysis  methodologies). 
Broadly  speaking,  implicit  or  explicit  procedures  are  available. 
Explicit  procedures  are  the  most  computationally  efficient  procedures 
since  they  do  not  require  (for  a  diagonal  mass  matrix)  equation  solving 
to  advance  the  solution.  However,  stability  restricts  the  size  of  the 
allowable  time  step.  On  the  other  hand,  unconditional  stability  can 
usually  be  achieved  in  implicit  procedures  but  they  do  require  solution 
of  a  system  of  equations  at  each  time  step.  First  and  foremost,  it  must 
be  pointed  out  that  a  purely  explicit  procedure  is  not  appropriate  for 
the  problem  at  hand  because  of  the  unreasonably  stringent  time  step 
restriction  resulting  from  the  presence  of  the  overly  stiff  fluid  in  the 
mixture  (even  for  highly  nonlinear  solid  material  models).  Methods  that 
combine  the  attractive  features  of  explicit  and  implicit  integrations 
have  recently  been  developed.  The  methods  used  here  fall  under  the 
category  of  "split  operator  methods."  In  operator  splitting  methods,  an 
implicit  integrator  is  selected  as  the  starting  point  and  the  integrand 
(right-hand  side  in  Equation  19)  is  split  so  that  the  system  of  equa¬ 
tions  solved  is  reduced.  The  specific  choice  made  is  obviously  problem 
dependent  as  discussed  further. 

Symbolically,  the  discretized  equations  of  motion  (Equation  19)  can 


be  written  as. 
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mSSm 


M  a  +  C  v  +  N(d , v) 


(26) 


where  M,  C,  N,  a,  v,  and  d  are  defined  by  Equation  19  in  an  obvious 
manner.  Time  integration  is  performed  by  using  the  implicit-explicit 
algorithm  of  References  35,  18,  and  19,  which  consists  of  satisfying  the 
following  equations: 


M  a  ^  +  C1  v  ^  +  CE  v  ,  +  NJ(d  _,,v  „,) 
~  ~n+l  -  ~n+l  ~  ~n+l  ~  ~n+l’~n+l 


.  i.E/j  ~  %  ..ext 

*  S  (VrVi>  *  ~n+i 


d  ,  =  d  .  +  o  At  a 

~n+l  ~n+l  ~n+l 


=  v  A  +  p  At  a 
~n+l  p  ~n+l 


where: 


=  d  +  At  v  +  (1  -  2p) 
~n  ~n 


=  v  +  (1  -  a)  At  a 
~n  ~n 


and  the  superscript  I  and  E  refer  to  the  parts  of  C  and  at  are  treated 


as  implicit  or  explicit,  respectively.  The  notation  is:  At  =  time 
ext  6xt 

step;  F  =  F  (t^);  and  an  are  the  approximations  to  d(tn), 

v(tn),  and  a(  )  I  a  and  p  =  algorithmic  parameters  that  control 
accuracy  and  stability  of  the  method.  It  may  be  recognized  that 


Equations  27  and  28  correspond  to  the  Newmark  formulas  (Ref  36).  The 


quantities,  dn+1  and  vn+^,  are  referred  to  as  "predictor"  values, 


while  d  and  v  are  referred  to  as  "corrector"  values.  From 
~n+l  ~n+l 


Equations  26  through  28,  it  is  apparent  that  the  calculations  are 


rendered  partly  explicit  by  evaluating  part  of  the  viscous  contri¬ 
bution  CE  v  and  the  force,  NE,  in  terms  of  data  known  from  the 

previous  step. 


Calculations  commence  with  the  given  initial  data  (i.e.,  dg  and 


Vq)  and  a^,  which  is  defined  by. 
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M  a0  =  Fjfc  -  CQ  vQ  -  N(d0,vQ) 


since  M  is  diagonal,  the  solution  of  Equation  28  is  rendered  trivial 


Impl ementati on 


At  each  time  step,  Equations  26  through  28  constitute  a  nonlinear 
algebraic  problem  that  is  solved  by  a  Newton-Ralphson  type  iterative 
procedure.  The  most  useful  and  versatile  implementation  is  to  form  an 
"effective  static  problem"  from  Equations  26  through  28  in  terms  of  the 


unknown  a  , ,  which  is  in  turn  linearized.  Within  each  time  step,  the 

~n+1  I  I 

calculations  are  performed  as  summarized  in  Table  1  in  which  C  and  K 


denote  the  parts  of  the  damping  and  tangent  stiffness  operators,  respec¬ 
tively,  to  be  treated  implicitly. 

The  following  choices  have  been  found  most  appropriate: 


.1 


Incompressible  Fluid  (Penalty  Formulation).  In  this  case,  C  is 
selected  to  contain  both  the  momentum  transfer  and  the  penalty  term 


2  +  C 


ss 
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(29) 
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where*  C 

~  (a,  p  =  s,  w)  =  damping  matrices  arising  from  the  penalty 
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treatment  of  the  fluid  contribution  as. 

P 


dQ 


(30) 


O' 


The  convective  fluid  force  (see  Equation  24)  is  usually  small  and  is 


treated  explicitly  with  no  resulting  computational  difficulty.  Note 


SW  W3  T  I 

that  since  C  =  (C  )  the  resulting  C  is  symmetric  for  the  choice 


adopted  here. 

As  for  the  solid  stress  force  (Equation  23)  contribution  to  the 
equations  of  motion,  three  options  are  possible:  "implicit,"  "explicit," 
or  "implicit-explicit"  treatment.  The  choice  is  made  as  follows: 
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Wave  Propagation  Type  Calculations .  Very  short  time  scale  (and 
high  frequency)  solutions  are  sought  and  an  explicit  treatment  of  the 
solid  effective  stress  contribution  is  usually  found  most  appropriate  in 
this  case.  The  time  step  size  restriction  resulting  from  stability 
considerations,  is  of  the  same  order  as  the  one  resulting  from  accuracy 
considerations  for  nonlinear  material  models.  Further,  the  specific 
implicit  treatment  adopted  for  the  fluid  contribution  allows  the  calcu¬ 
lations  to  be  carried  out  at  a  time  step  usually  close  to  the  time  step 
corresponding  to  the  propagation  of  the  solid  compressional  wave  through 
the  solid  phase  of  the  critical  element. 

Vibration  Type  Calculations.  Since  the  frequencies  captured  are 
usually  much  higher  than  above,  an  implicit  treatment  of  the  solid 
effective  stress  contribution  is  usually  convenient  in  this  case  since 
it  allows  the  time  step  to  be  selected  following  accuracy  considerations 
only.  Unconditional  stability  is  achieved  by  selecting  the  proper 
algorithmic  parameters  as  discussed  later.  However,  for  nonlinear 
analyses,  a  purely  implicit  treatment  requires  a  maxtrix  reform/ 
factorize  at  each  time  step  (and  for  every  iteration  to  be  performed,  in 
general),  thus  producing  a  considerable  computational  burden.  It  is 
therefore  convenient,  in  the  nonlinear  case,  to  adopt  an  implicit- 
explicit  treatment  of  the  effective  stress  contribution  as  follows:  The 
linear  part  of  the  stiffness  is  treated  implicitly  while  the  remaining 
nonlinear  part  is  treated  explicitly.  For  that  purpose,  a  solid  stiff¬ 
ness  operator  is  defined  from  Equations  9  and  23  through  linearization 
as : 


where:  KS  =  material  tangent  part 

G 

K  =  "initial„stress"  or  geometric  part,  formed  from  the  tensors 
DS  and  Db  (Equation  9)  in  the  usual  manner  (see  Ref  16). 


Table  1.  Flow  Chart 


Ini tial i zation: 
i  =  0 


Q1  =  Q2  =  0 

Predictor  Phase: 
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Form  Effective  Mass:  (Reform  and  factorize  only  if  required) 


M*  =  M  +  At  a  C1  +  At2  p  K1 


Form  Residual  Force: 
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In  the  implicit-explicit  procedure: 


(32) 


where:  Ks  =  linear  elastic  contribution  to  the  material  tangent 

stiffness  (from  Equation  10). 

Such  a  choice  does  not  always  lead  to  unconditional  stability.  The 

G 

difficulty  is  not  usually  associated  with  the  explicit  treatment  of  K 
(which  contains  terms  of  the  stress  order  and  therefore  usually  has  a 
negligible  impact  on  stability),  but  rather  from  the  explicit  treatment 
of  the  nonlinear  term  (KS  K ^)  for  materials  with  a  locking  tendency. 
In  that  case,  care  must  be  exercised  in  selecting  a  time  step  smaller 
than  the  one  associated  with  the  fastest  expected  wave  speed  correspond¬ 
ing  to  the  subsequent  stress  histories  to  be  followed  by  the  material 
e 1 ement  s 


Diffusion  Type  Calculations .  It  is  sometimes  desirable  to  capture 
the  purely  diffusive  part  ("consolidation"  part)  of  the  solution 
dynamically  "  Such  a  necessity  arises  in  situations  in  which  both  short 
and  1  onq  time  solutions  to  a  dynamical  problem  are  sought  (such  as  in 
seismic  or  blast  induced  liquefaction  simulations).  As  shown  in  Refer¬ 
ences  14  and  16,  by  switching  to  an  appropriate  choice  of  the  Newmark 
parameters,  n  =  3/2  and  ft  =  1 ,  and  by  using  the  implicit-explicit  option 
described  above,  all  dynamic  transients  can  be  damped  out,  and  purely 
diffusive  ( conso 1 i dat ’ on )  solutions  can  be  obtained  "dynamically"  by 
solving  the  dynamic  equations 

The  penalty  treatment  of  the  storage  equation  constraint  requires 
that  \w  be  selected  as  a  "large  number"  capable  of  predominating  the 
It  should  be  picked  according  to  the  relation: 


other  modu  '  i 


where:  C  =  a  constant  that  depends  only  on  the  computer  word  length. 
Numerical  studies  reveal  that  for  floating-point  word  lengths  of  60  to 
64  bits,  C  ~  107. 

It  has  been  determined,  on  the  basis  of  numerical  experiments,  that 
it  is  better  to  use  reduced  integration  of  the  penalty  terms  and  of  the 
solid  volumetric  stiffness  (div  vS  contribution)  to  alleviate  potential 
mesh  locking  phenomena. 

Compressible  Fluid.  In  this  case,  C*  contains  the  momentum  trans¬ 
fer  term  contribution  to  the  equations  of  motion,  viz. , 


-Z  -z 


and  the  convective  fluid  force  is  treated  explicitly. 

Again,  the  fluid  pressure  contribution  is  treated  implicitly.  For 
this  purpose,  a  fluid  stiffness  operator  is  defined  from  Equations  14 
and  19  through  linearization  as 


Where  C  p  (a,  p  =  s,  w)  are  the  same  matrices  as  previously  defined  in 
Equation  31.  However,  note  that  in  this  case,  the  matrices  contribute 
to  the  stiffness  matrix  rather  than  to  the  damping  matrix. 

The  same  options  as  described  previously  for  the  solid  stress  force 
contribution  to  the  equations  are  used.  Note  that  in  the  implicit- 
explicit  procedures, 


K*  +  CSS  CSW 


which  again  is  a  symmetric  matrix. 
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Stabi 1 i ty.  The  resulting  stability  conditions  (Ref  35)  are  sum¬ 
marized:  In  all  cases  a  >  1/2: 

1.  Implicit  Treatment:  Unconditional  stability  is  achieved  if  p  > 
or/2  and  it  is  recommended  that  (Ref  37): 


(37a) 


Implicit-Explicit  Treatment:  The  time  step  restriction  is 


w  At  < 


(37b) 


to  maximize  high-frequency  numerical  dissipation.  ,  w  =  highest  natural 
frequency  associated  with  the  explicit  part  of  the  stiffness  operator. 
The  maximum  expected  frequency  may  be  bounded  by  the  frequency  of  the 
smallest  element,  viz.,  for  a  rectangular  four  node  bilinear  element, 


(37c) 


where:  L  =  smallest  dimension  of  the  element 

C  =  wave  speed  (  =  \(AS  +  2  pS)/pS  for  the  linear  model). 


Program  DYNAFLOW 


Equations  19  through  37  have  been  incorporated  into  the  finite 
element  computer  program  DYNAFLOW  (Ref  38).  Numerical  results  that 
illustrate  the  performance  of  the  proposed  numerical  schemes  in 
analyzing  the  propagation  of  plane  progressive  waves  (Ref  39,  1,  29,  26, 
40,  and  41)  in  f 1 uidsaturated  porous  soil  media  have  been  reported  in 
Reference  13.  In  Reference  13  both  compressible  and  incompressible 
fluid  cases  are  considered.  In  the  following,  attention  is  restricted 
to  vibration  calculations  in  soil-structure  interacting  systems  associa¬ 
ted  with  seismic  events  in  which  the  pore  fluid  may  be  assumed  to  be 
incompressible  ( compared  to  the  compre.-sibi  1  ity  of  the  soil  skeleton). 


In  the  calculation  reported  the  four-node  bilinear  isoparametric 

element  (see  Ref  43  and  33  for  a  detailed  description)  is  used  with  the 

standard  selective  integrations  scheme  (Ref  43).  Also,  the  Newark 

algorithmic  parameters  are  always  selected  such  that  a  >  1/2  and 

2 

p  =  (a  +  1/2)  /4  to  maximize  high-frequency  numerical  dissipation 
(Ref  37). 

Time  integration  of  the  semi -discrete  finite  element  equations  is 
performed  by  using  the  implicit-explicit,  predictor-(multi )  corrector 
option.  The  structural  domains  and  the  fluid  phase  are  always  treated 
implicitly,  whereas  the  solid  skeleton  phase  is  treated  partly 
implicitly,  partly  explicitly.  Specifically,  the  "elastic"  part  of  the 
solid  stiffness  is  treated  implicitly,  whereas  the  elastic-plastic  part 
of  the  stiffness  is  treated  explicitly.  Such  a  choice  allows  time 
integration  to  be  performed  most  efficiently  for  cases  in  which  many 
load/unload  cycles  occur,  such  as  during  dynamic  events. 

In  order  to  simulate  realistic  initial  conditions  for  the  stresses 
and  strains  in  the  soil  deposits  before  input  of  the  earthquake  like 
ground  motions,  in  all  the  example  problems  reported,  the  numerical 
simulation  was  carried  in  at  least  two  sequential  steps: 

1.  Installation  of  Soil  Deposit.  This  step  in  the  numerical 
simulation  was  designed  to  simulate  realistic  initial  conditions.  The 
computed  gravitational  stresses/displacements/pore-water  pressures  were 
used  as  initial  conditions  for  the  ground  motion  calculations. 

2.  Ground  Motion  Simulation:  The  computed  gravitational  stresses/ 
displacements/pore-water  pressures  computed  in  the  first  sequence  were 
used  as  initial  conditions  for  the  ground  motion  calculations.  In  order 
to  avoid  the  initial  propagations  of  any  spurious  noises,  the  accelera¬ 
tion  array  was  first  cleared  before  activating  the  ground  motion.  Also, 
the  algorithmic  parameters  were  reset  to  a  =  0.65,  p  =  0.33  (slightly 
diffusive),  the  time  to  zero,  and  the  time  step.  At,  selected  to 
properly  follow  the  details  of  the  stresses  (and  pore-water  pressures) 
"at  rest"  before  input  of  the  ground  shaking.  For  that  purpose,  small 
initial  stresses  were  first  input  in  the  soil  deposit  and  gravity  was 


turned  on.  A  number  of  steps  (and  iterations)  were  used  to  apply 
gravity  with  a  large  enough  time  step.  At,  to  ensure:  (1)  that  no 
excess  pore-water  pressures  (over  steady  state  conditions)  would  build¬ 
up.  ,  (2)  full  development  of  effective  stresses  in  the  foundation  soil 
would  occur;  and  (3)  that  no  transients  would  occur.  For  that  purpose, 
a  backward  scheme  with  high  numerical  dissipation  was  used  to  damp  out 
all  transients  setting  a  =  1.50  and  0  =  1.00.  Iterations  (Newton- 
Raphson  or  modified  Newton-Raphson)  were  used  to  ensure  proper  con¬ 
vergence  of  the  solution  at  each  time  step. 

Post-Dynamic  Event  Simulation 

In  saturated  soil  deposits,  excess  pore-water  pressures  usually 
buildup  and  do  not  have  time  to  fully  dissipate  within  the  time  frame  of 
the  earthquake.  Subsequent  redistribution  and  diffusion  of  these  excess 
pore-water  pressures  following  the  dynamic  event  are  often  of  interest 
since  they  may  lead  to  failures  (e.g. ,  in  earth  dams  where  failures  have 
been  recorded  in  the  field  several  minutes  to  several  hours  following 
the  earthquake).  It  is,  however,  usually  inefficient  (costwise)  to 
capture  numerically  the  post-event  behavior  with  a  time  step  designed  to 
follow  the  details  of  an  earthquake  motion.  Therefore,  in  the  follow¬ 
ing,  post-dynamic  events  have  been  calculated  by  resetting  a  =  1.5  and 
0  =  1.00  (highly  diffusive  backward)  and  by  selecting  the  appropriate 
time  step,  At,  for  capturing  the  post-event  diffusion  part  of  the 
solution.  (Note.  Either  constant  At  or  geometrically  growing  At  were 
used  for  this  sequence.) 

PRINCETON  UNIVERSITY  EFFECTIVE  STRESS  SOIL  MODEL 

Theory 


The  hysteretic  stress-strain  behavior  of  the  soil  skeleton  modeled 
by  using  the  effective-stress  elastic-plastic  model  reported  by 
Professor  Prevost  (Ref  31).  The  model  is  an  extension  of  the  simple 


multi-surface  ^-plasticity  theory  and  uses  conical  yield  surfaces.  The 
model  has  been  developed  for  cohesionless  soils  and  tailored  to: 

1.  Retain  the  extreme  versatility  and  accuracy  of  the  multi¬ 
surface  J.-theory  in  describing  observed  shear  nonlinear 
hysteretic  behavior,  shear  stress- induced  anisotropy  effects. 

2.  Reflect  the  strong  dependency  of  the  shear  dilatancy  on 
the  effective  stress  ratio  in  granular  cohesionless  soils. 

The  model  is  applicable  to  general  three-dimensional  stress-strain 
conditions,  but  its  parameters  can  be  derived  entirely  from  the  results 
of  conventional  triaxial  soil  tests. 

Yield  Function.  The  yield  function  selected  is  of  the  following 

form: 

f  =  |(s  -  p  a):  (s  -  p  or)  -  m2  p2  =  0  (38) 

where:  s  =  ct  -  p6  =  deviatoric  stress  tensor 

<v  /V 

p  =  1/3  tr  6  =  effective  mean  normal  stress 
6  =  effective  stress  tensor 

a  =  kinematic  deviatoric  tensor  defining  the  coordinates 
of  the  yield  surface  center  in  deviatoric  stress 
subspace  ("the  back  stress") 

m  =  material  parameter 

tr  =  trace 

The  yield  function  plots  as  a  conical  yield  surface  (Drucker-Prager)  in 
stress  space,  with  its  apex  at  the  origin,  and  is  shown  in  Figure  1. 

Unless  a  =  0,  the  axis  of  the  cone  does  not  coincide  with  the  space 

diagonal.  The  cross  section  of  the  yield  surface  by  any  deviatoric 
plane  (p  =  constant)  is  circular.  Its  center  does  not  generally  coin¬ 
cide  with  the  origin,  but  is  shifted  by  the  amount  p  a.  This  is 
illustrated  in  Figure  1  in  the  principal  stress  space. 
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Figure  1.  Yield  surfaces  in  principal  stress  space. 


Flow  Rule.  The  plastic  potential  is  selected  such  that  the 
deviatoric  plastic  flow  be  associative.  However,  a  non-associative  flow 
rule  is  used  for  its  dilatational  component,  and  in  the  following: 

D  =  Inn  l 


where:  tr  =  trace 

P  =  symmetric  second-order  tensor  which  defines  (in  stress- 
space)  the  "direction"  of  plastic  deformations 

q  =  (3/2  s:s)  /p  =  stress  ratio 

H  =  material  parameter. 

When  n<H,  the  trace  of  P  is  less  than  zero  and  plastic  compaction  takes 

place,  whereas  when  n<rl*  the  trace  of  P  is  greater  than  zero  and  plastic 

dilation  takes  place.  The  case  n  =  H  corresponds  to  no  plastic 

-  -  3  - 

volumetric  strains.  In  the  following,  n  =  when  tr  s  <0,  and  n  =  Or 

3  L  ~  t 

when  tr  s  >0. 


Hardening  Rules.  In  order  to  allow  for  the  adjustment  of  the 
plastic  hardening  rule  to  any  kind  of  experimental  data,  a  collection  of 
nested  yield  surfaces  is  used.  The  yield  surfaces  are  all  similar 
conical  surfaces  (Equation  1),  and  a  purely  kinematic  hardening  rule  is 
adopted.  Upon  contact,  the  yield  surfaces  are  translated  by  the  stress 
point  and  the  direction  of  translation  is  selected  so  that  no  over¬ 
lapping  of  the  surfaces  can  take  place  (see  Reference  31  for  more 
details). 

Each  yield  surface  is  associated  with  a  plastic  modulus  H'.  The 
dependence  of  the  plastic  modulus  upon  the  deviatoric  stress  is  assumed 
to  be  in  the  following  form: 


H' 

where: 


H'  -  H' 
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H'  +  H1 
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s  =  3(s  -pa) 


(40) 


and  =  material  parameters. 

Effects  of  the  Effective  Mean  Confining  Stress.  The  dependence  of 
the  material  moduli  upon  the  effective  mean  confining  stress  is: 


G  = 


=  HI 


(41) 


where: 


G  =  elastic  shear  modulus 


B  =  elastic  bulk  modulus 
H'  =  plastic  modulus 

n  =  experimental  parameter  (n  =  0.5  for  most 
cohesionless  soils) 


r 


=  moduli  at  reference  effective  mean  normal  stress  p 


LABORATORY  TRIAXIAL  TEST  SIMULATION  -  COHESIONLESS  SOILS 


Test  Program 

To  validate  the  Princeton  University  Effective  Stress  Soil  Model,  a 
series  of  sands  were  evaluated  in  both  drained  and  undrained  triaxial 
tests.  The  tests  were  carried  out  in  both  extension  and  compression, 
and  over  a  range  of  confining  cell  pressures  and  relative  densities.  To 
eliminate  some  of  the  uncertainty  as  to  the  effects  of  sample  prepara¬ 
tion  and  testing  procedure  on  the  model's  sensitivity,  the  test  data 
were  obtained  from  different  university  experts  and  the  Waterways 
Experiment  Station. 

The  material  parameters  for  the  Princeton  University  soil  model  are 
derived  from  drained  triaxial  compression  and  extension  data,  defined  by 
shear  stress,  shear  strain,  mean  stress,  and  volume  strain  information. 
This  information  forms  the  input  for  the  computer  Program  MUD,  which  is 
an  automated  procedure  to  derive  the  material  model  moduli  and  yield 
surface  parameters.  The  material  model  has  been  implemented  into  a 
one-point  equilibrium  solution  procedure,  Program  TEST  A,  to  provide  a 
quick  and  easy  means  to  evaluate  a  particular  soil  model  in  various 
stress  states  for  direct  comparisons  and  evaluation  of  laboratory  tests. 
TEST  A  was  used  for  this  portion  of  the  evaluation,  which  is  the  simula¬ 
tion  of  the  triaxial  test  data. 

The  basic  steps  performed  to  evaluate  the  laboratory  test  data 

were : 


1.  Select  a  representative  set  of  drained  triaxial  tests  for  which 
both  compression  and  extension  results  are  available. 

2.  Use  the  results  from  step  one  as  input  to  Program  MUD  to  deter¬ 
mine  model  parameters. 

3  Construct  the  TEST  A  input,  with  the  output  from  Program  MUD  to 
evaluate  the  models  representation  of  input  test  data. 
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4.  Analyze  the  stress  ratio  value  in  the  model  parameters  to 
determine  a  reasonable  estimate  of  the  critical  stress  ratio  value.  The 
critical  stress  ratio  parameters  control  the  compactive-dilative  char¬ 


acteristics  of  the  model,  and  thus  the  volumetric  strain. 

5.  Adjust  the  critical  stress  ratios  through  successive  TEST  A 
analysis  to  optimize  the  drained  model  versus  undrained  prediction 
agreement. 

6.  Once  soil  model  agreement  is  achieved  with  the  original  test 
data,  the  parameters  determined  to  quantify  the  soil  behavior  may  be 
used  to  predict  behavior  at  different  confining  pressures  and  along 
different  stress  paths  without  further  modification. 

Monterey  "0"  Sand 

The  University  of  California,  Davis,  under  contract  with  the  Naval 
Civil  Engineering  Laboratory  (NCEL)  performed  triaxial  tests  on  Monterey 
"0"  sand.  Figures  2(a)  and  2(b)  present  the  drained  compression  and 
extension  results  from  the  triaxial  test  data  that  were  used  to  formu¬ 
late  the  model,  and  the  generated  results  from  Program  TEST  A.  As 

noted,  agreement  is  nearly  exact.  In  Figures  3(a)  and  3(b)  the  model 

2 

parameters  determined  at  an  initial  confining  stress  of  1  kg/cm  are 

used  to  simulate  drained  compression  and  extension  tests  at  an  initial 

2 

confining  stress  of  3  kg/cm  .  Again,  the  agreement  between  predicted 
versus  observed  behavior  is  very  good. 

Figures  4(a)  and  4(h)  present  the  shear  stress-strain  and  effec¬ 
tive  stress  path  for  an  undrained  compression  test.  The  undrained 
stress  path  agrees  almost  perfectly,  which  is  a  direct  indication  of  the 
model's  ability  to  capture  correctly,  the  soils  compactive-dilative 
characteristics.  This  in  turn  generates  the  positive  and  negative  pore 
pressures,  and  hence  the  changes  in  effective  stress,  which  agree  very 
well  with  those  observed  in  the  actual  test.  In  Figure  4(a)  the  shear 
stress-shear  strain  comparisons  for  the  undrained  analysis  are  again  in 
good  agreement. 


computed 

observed 


s 

/ 


/ 

/ 

f 


0.^0  O.dO  0.C0  0.30  1.00  l.icQ 

Shear  Strain  (cm/cm) 

(a)  Compression 


i  X  !  0  1 1 

i  -  -  -  i 

\ .  >io  i .  r,o 


I 


i 


- -  computed 

- observed 


f  f-  i  - 

-J  CO  J  i)u  c  .so 

Shear  Strain  tun  cm > 

(b>  lx  tension 


f 

t  .  30 


Figure  2.  Drained  compression  and  extension  Monterey  "0"  sand, 
confining  cell  pressure  =  1  kg/cm2. 
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Figures  5,  6(a)  and  6(b)  are  a  demonstration  of  the  same  material 
model  derived  from  drained  monotonic  tri axial  data,  predicting  an 
undrained  cyclic  triaxial  test.  As  can  be  seen  in  Figure  5,  the  model 
does  an  adequate  job  of  predicting  the  observed  hysteretic  behavior. 
There  is  a  minor  difference  in  initial  stiffnesses  between  predicted  and 
observed  behavior  caused  by  a  discrepancy  in  tracking  initial  volumetric 
test  data.  This  could  have  been  improved  by  providing  added  definition 
to  the  initial  portion  of  the  test  data  but  was  judged  satisfactory  as 
it  is.  Note  the  general  agreement  between  the  stiffnesses  in  both 
loading  and  unloading,  which  demonstrates  correct  energy  dissipation  in 
each  cycle.  In  the  stress  path  plot  in  Figure  6(a)  the  same  degradation 
per  cycle  is  demonstrated.  Note  the  number  of  predicted  cycles  to 
failure  is  4.5,  this  agrees  well  with  the  observed  5.5  cycles.  The 
overall  generation  of  pore  water  pressure  is  in  good  agreement  as  well, 
which  is  shown  both  in  the  decline  in  effective  stress  to  failure  in 
Figure  6(a),  as  well  as  in  the  maximum  level  attained,  Figure  6(b). 
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Figure  5.  Undrained  cyclic  analysis  Monterey  "0"  sand  confining 
cell  pressure  =  1  kg/cm2,  shear  stress-vertical  strain 
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Figure  6.  Undrained  cyclic  analysis  Monterey  "0"  sand,  confining 
cell  pressure  =  1  kg/cm2,  stress  path  and  pore  pressure. 


Leighton-Buzzard  Sand 


Triaxial  tests  on  Leighton-Buzzard  sand  were  performed  by  the 
University  of  California,  Davis,  under  contract  to  NCEL.  Figures  7(a) 

and  7(b)  show  the  drained,  monotonic  compression  and  extension  com¬ 
parisons.  Figure  8(a)  shows  the  undrained  shear  stress-shear  strain 

data,  and  Figure  8(b)  shows  the  undrained  stress  path  from  the  same 

test.  The  model  in  this  example  predicts  a  softer  initial  response  with 
a  smaller  generation  of  pore  water  pressure  than  the  recorded  test  data, 
but  these  differences  are  well  within  the  ranges  of  acceptable  differ¬ 
ences  between  two  laboratory  tests. 

The  simulation  of  the  undrained  cyclic  test  from  the  monotonic 
model  again  shows  good  agreement  in  producing  an  accurate  assessment  of 
energy  dissipation  and  failure  state.  The  model  in  this  case  predicted 
1.75  cycles  to  failure  compared  to  the  observed  failure  at  2.75, 
(Figures  9,  10(a)  and  10(b)).  The  final  effective  stress  values  and 

pore  water  pressures  are  in  good  agreement  as  well. 

Silica  Sand 

Silica  sand  was  tested  at  the  University  of  California  at  Los 

Angeles,  under  contract  with  NCEL.  Figures  11(a)  and  11(b)  show  the 

reproduction  of  the  drained  triaxial  compression  and  extension  data  used 

to  generate  the  material  model.  As  would  be  expected  agreement  is  very 

good.  To  demonstrate  the  model's  ability  to  predict  tests  at  different 

2  2 

confining  cell  pressures,  simulations  were  run  at  5  kg/cm  and  1  kg/cm  , 
the  results  are  plotted  on  Figures  12(a)  and  12(b),  respectively.  Again 
generated  curves  are  in  very  good  agreement  with  the  observed  results. 

Figures  13(a)  amd  13(b)  show  the  model's  ability  to  simulate  a  sub¬ 
stantially  different  stress  path.  The  data  plotted  in  these  figures  is 
from  a  drained  uniaxial  test  (kQ)  where  vertical  stress  is  increased 
while  holding  horizontal  strain  at  zero.  The  model  uses  the  same  para¬ 
meters  used  to  predict  the  standard  drained  and  undrained  tests  for  the 
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silica  sand.  Shear  stress-vertical  strain  agreeeient  and  volume  stress- 
volume  strain  agreement  are  very  good.  In  addition  to  the  uniaxial 
test,  a  proportional  loading  test  was  run.  The  results  are  presented  in 
Figures  14(a)  and  14(b).  The  initial  computed  stress  path  deviates 
slightly  from  the  observed  data,  however,  the  results  do  present  a 
reasonable  simulation  of  the  actual  test  once  the  initial  anisotropic 
stability  is  reached. 

Figures  15,  16(a),  and  16(b)  present  an  undrained  cyclic  triaxial 
simulation  using  the  same  drained  monotonic  model  used  throughout  the 
analysis  of  the  silica  sand.  The  model  predicts  liquifaction  in  4.7 
cycles  as  compared  to  measured  values  of  5.7  and  3.4  in  two  identical 
laboratory  tests.  Rates  of  energy  dissipation  and  pore  water  pressure 
generation  are  again  in  reasonable  agreement,  and  as  in  prior  analysis, 
the  failure  state  of  the  soil  is  captured  quite  well. 

In  addition  to  the  determination  of  model  parameters  from  the 
consolidated  drained  tests  where  volume  strains  are  measured,  analyses 
were  performed  to  predict  drained  and  undrained  tests  from  an  undrained 
model.  In  an  undrained  test,  the  volume  strain  is  held  constant  (i.e., 
to  zero)  and  excess  pore  water  pressure  is  developed.  This  pore  water 
pressure  is  a  result  of  the  skeletal  deformation,  which  is  an  indication 
of  the  degree  of  skeletal  volume  strain. 

Program  MUD  was  used  with  undrained  test  data  to  develop  the  model 
parameters.  (Previously  drained  data  was  used.)  These  parameters  were 
then  used  to  predict  drained  tests.  The  results  of  the  undrained  model 
for  silica  sand  are  presented  in  Figures  17(a)  and  17(b).  The  compres¬ 
sion  and  extension  results  both  agree  very  well  with  the  actual  test 
data,  and  the  previous  results  calculated  using  the  model  formulated 
from  drained  data  (Figures  11(a)  and  11(b)). 

The  undrained  model  predicting  undrained  compression  is  in  good 
agreement  with  the  data  and  drained  model  predictions.  There  is  a 
slightly  larger  pore  pressure  build  up  in  the  undrained  model  (Fig¬ 
ures  17(a)  and  17(b))  shown  by  the  larger  decrease  in  effective  stress 
presented  in  the  stress  path  plot  (Figure  17(b)). 
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Figure  8.  Undrained  compression  Leighton-Buzzard  sand  confining  cell 
pressure  =  1  kg/cm2. 
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Figure  10.  Undrained  cyclic  analysis  Leighton-Buzzard  sand,  confining 
cell  pressure  =  1.0kg/cm2,  stress  path  and  pore  pressure. 
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Figure  11.  Drained  compression  and  extension  silica  sand,  confining 
cell  pressure  =  2  kg/cm2. 
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Figure  14.  Proportional  loading  test  with  silica  sand. 
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Figure  15.  Undrained  compression  silica  sand,  confining  cell 
pressure  =2.0  kg/cm2. 
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Figure  16.  Undrained  cyclic  analysis  silica  sand,  confining  cell 
pressure  =2.0  kg/cm2,  shear  stress-verticcl  strain. 
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Figure  17. 


Undrained  cyclic  analysis,  silica  sand,  confining  cell 
pressure  =  2.0kg/cffl2,  stress  path  and  pore  pressure. 
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Miesers  Bluff  Sand 


The  Miesers  Bluff  sand  was  tested  at  the  Waterways  Experiment 
Station  Geotechnical  Laboratory.  Figures  18(a)  and  18(b)  show  the 
drained  test  data  used  to  generate  the  model  and  the  computed  results. 
Figures  19(a),  19(b),  and  20  present  the  model's  ability  to  predict  the 
drained  response  of  the  material  over  a  range  of  confining  pressures 
from  0.14  MPA  to  6.9  MPA  based  upon  the  model  parameters  determined  at 
3.45  MPA.  Again  agreement  is  very  good  across  the  entire  range  of  the 
drained  analysis.  Unfortunately,  many  of  the  reference  curves  lack 
definition,  because  the  test  data  used  for  this  series  had  to  be 
digitized  from  hard  copy  plots. 

Figures  21(a)  and  21(b)  present  an  undrained  simulation  run  at  a 
confining  stress  of  3.45  MPA.  Note  that  the  trends  of  this  particular 
sand  are  well  modeled  with  reasonable  accuracy,  but  due  again  to  the 
coarseness  of  the  digitizing  process,  both  the  test  data  and,  therefore, 
model  predictions  of  stress  strain  behavior  may  be  somewhat  in  error. 
Figures  22(a)  and  22(b)  present  the  results  of  an  additional  undrained 
test  run  at  1.72  MPA.  Again,  note  the  agreement  between  the  predicted 
and  observed  behavior. 

Figures  23  and  24  represent  an  undrained  cyclic  triaxial  simulation 
of  the  Miesers  Bluff  sand.  Cyclic  test  data  were  not  available  for  this 
analysis  but  the  simulation  was  included  for  completeness. 


SENSITIVITY  OF  PARAMETERS 

Critical  Stress  Ratio 

The  critical  stress  ratio  (CSR)  is  the  parameter  in  extension  or 
compression  which  determines  the  stress  state  at  which  there  is  a  change 
between  compactive  and  dilative  volume  strains.  Adjustment  of  this 
parameter  has  little  effect  on  the  shear  s tres s- s t ra i n  curves  in  the 
drained  simulation,  but  does  make  significant  changes  in  the  amount  of 
volume  strain  generated  by  the  model  An  example  of  varying  the  critical 


stress  ratio  by  ±6%  is  presented  for  a  drained  test  using  the  silica 
sand  model  (Figures  25(a)  and  25(b)). 

Figures  26(a)  and  26(b)  present  the  results  of  the  same  modifi¬ 
cation  in  the  critical  stress  ratio  during  an  undrained  simulation. 
Note  the  significant  changes  in  strength  caused  by  the  10%  increase  in 

critical  stress  ratio  leading  to  a  prediction  of  near  total  loss  of 

strength.  Decreasing  the  value  of  the  stress  ratio  by  6%  more  than 

doubles  the  amount  of  shear  strength  developed.  Figure  26(b)  presents 
the  stress  paths  associated  with  each  of  the  modifications.  Note  that 
each  of  the  models  follows  essentially  the  same  path  up  to  the  point  of 
contacting  the  critical  state  line.  At  this  point,  if  the  sample  has 
not  begun  to  dilate,  it  then  loses  strength  by  moving  toward  zero 

effective  stress. 


Elastic  Moduli 

Figures  27(a)  and  27(b)  demonstrate  the  effects  of  varying  the 
elastic  moduli  for  all  surfaces  in  the  drained  material  model  by  ±25%. 
The  effects  on  the  shear  stress-strain  curves  are  minimal;  however, 
there  is  a  more  significant  change  in  the  degree  of  compactive  volume 
strain  generated.  Note  as  well,  the  changes  in  volume  strain  are  not 
equivalent  for  both  cases.  This  is  caused  by  an  increase  in  plasticity 
associated  with  the  decreased  elastic  stiffness  necessary  to  reach  the 
same  load  level.  The  increased  elastic  stiffness  in  turn,  allows  the 
given  load  to  be  reached  at  a  smaller  strain  level,  but  does  not 
decrease  the  plastic  strain  to  the  same  degree.  Figures  28(a)  and  28(b) 
present  the  effects  of  the  same  variation  on  the  undrained  simulation. 

Plastic  Moduli 


Varying  the  plastic  moduli  for  all  surfaces  has  significant  effects 
on  both  the  drained  and  undrained  simulations  (Figures  29  and  30).  A 
greater  change  appears  as  a  result  of  increasing  the  plastic  moduli  for 
each  yield  surface  than  by  decreasing  it.  The  effect  is  two-fold  in 


that  it  not  only  increases  the  overall  stiffness  of  the  system  in  the 
drained  simulation  (Figure  29(a)),  but  decreases  the  overall  compactive 
volumetric  strain  as  well  (Figure  29(b)).  These  factors  both  contribute 
to  the  greatly  enhanced  stiffness  in  the  undrained  simulation 
(Figure  30(a))  by  improving  the  skeletal  response  while  at  the  same  time 
decreasing  the  excess  pore  water  pressure  generation. 

Reducing  the  plastic  moduli  has  less  et feet  in  both  drained  and 
undrained  simulations  causing  a  slightly  larger  amount  of  compaction 
before  dilation  begins.  This  produces  a  moderately  softer  system  in 
both  the  drained  and  undrained  environments. 

Yield  Surface  Size 

Changing  the  yield  surface  size  changes  the  shear  strength  values 
which  may  be  achieved  at  a  particular  strain  level.  Figures  31  and  32 
clearly  demonstrate  this  property  of  the  model. 

In  the  drained  analysis  (Figures  31(a)  and  31(b)),  increasing  the 
yield  surface  size  decreases  the  shear  strain  at  a  particular  level  of 
shear  stress,  and  reduces  volumetric  strains  not  by  a  change  in  moduli 
values,  but  rather  as  a  result  of  changing  the  position  of  compactive- 
dilative  interchange  in  relation  to  the  original  critical  surface 
position.  Reducing  the  yield  surface  size  results  in  the  reaching  of 
the  various  surfaces  earlier,  at  reduced  stress  levels,  and  induces 
greater  amounts  of  plasticity.  At  the  same  time,  the  reduction  also 
moves  the  compacti vedi 1 ati ve  point  outside  the  outermost  surface, 
allowing  large  compactive  volumetric  strains  to  occur  (Figure  31(b)). 

In  the  undrained  simulation  (Figures  32(a)  and  32(b))  the  two-fold 
effect  of  a  stiffer  skeleton  and  reduced  volumetric  strains  are  present 
The  increase  in  yield  surface  size  produces  a  stiffer  system  with  a 
reduction  in  positive  pore  water  pressure  generation  before  dilation 
begins.  The  reduction  in  size  causes  an  early  generation  of  volume 
strains  and  positive  pore  water  pressures,  coupled  with  a  softer  system 
resulting  in  the  zero  effective  stress  condition 
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Figure  26.  Critical  stress  ratio  variation  ±6%,  +10%,  undrained 
compression  shear  stress-strain,  stress  path. 
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Yield  Surface  Position 


Variation  of  yield  surface  position  has  equal  effects  in  both 
drained  and  undrained  simulations  (Figures  33  and  34).  By  moving  the 
surfaces  in  a  positive  direction,  the  system  is  stiffened  by  changing 
the  points  at  which  the  load  path  intercepts  -  particular  surface, 
i.e.,  the  path  achieves  greater  stress  levels  with  less  plasticity 
before  impacting  the  next  surface,  thereby,  remaining  at  a  given  stiff¬ 
ness  through  a  larger  strain  increment.  By  moving  the  surface  position 
in  a  negative  direction,  the  load  path  then  intercepts  a  particular 
surface  at  a  smaller  stress  increment,  and  induces  greater  amounts  of 
plasticity,  effectively  softening  the  system. 

Effects  on  the  undrained  simulation  are  presented  in  Figures  34(a) 
and  34(b)  Qualitatively,  the  changes  produce  the  same  response  as 
demonstrated  in  the  drained  simulation. 


LABORATORY  TESTS  -  COHESIVE  SOILS 

Kao  1 i ni te  I 


The  triaxial  tests  for  the  kaolinite  clay  fitted  in  Figure  35  were 
performed  at  the  University  of  California,  Davis.  The  model  parameters 
for  the  kaolinite  were  derived  from  the  consolidated  undrained  triaxial 
compression  and  extension  tests.  Figures  35(a)  and  35(c)  present  the 
predictions  for  shear  stress  versus  shear  strain.  Th^re  is  some  dis¬ 
crepancy  in  the  predicted  ultimate  strength  in  compression  which  may  be 
attributed  to  the  error  involved  in  digitizing  the  data  points  for  the 
effective  volume  stress  based  on  the  stress  path  plots  of  the  original 
data  The  extension  predictions  in  Figure  35(c)  are  in  excellent 
agreement  with  the  measured  values.  The  stress  path  predictions 
presented  in  Figures  35(c)  and  35(d)  show  the  same  behavior  as  were 
demonstrated  in  the  shear  stress  strain  plots,  with  acceptable  agreement 
tn  both  cases  Better  agreement  may  be  achieved  by  increasing  the 
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number  of  yield  surfaces  in  the  initial  portion  of  the  loading.  This 
requires  better  data  definition  and  a  much  finer  digitizing  process. 

Kaolinite  II 

The  soil  data  used  in  this  section  are  part  of  the  data  which  were 
collected  by  the  organizing  committee  of  the  National  Science  Foundation 
North  American  Workshop  in  soil  engineering  held  May  28-30,  1980  at 
McGill  University,  Montreal,  Canada  and  reported  in  Reference  40.  Axial 
test  data  on  the  laboratory  prepared  kaolinite  clay  had  been  provided. 
Predictions  about  the  constitutive  behavior  of  the  soil  subjected  to 
loading  stress  paths  not  identified  in  the  data  had  been  requested  and 
their  analysis  compares  the  model  predictions  with  observed  behavior  in 
the  tests. 

The  dashed  lines  in  Figure  36  show  the  experimental  results 
obtained  in  conventional  triaxial  undrained  monotonic  axial  compression/ 
extension  soil  tests  and  the  solid  lines  show  the  design  curves  used  to 
determine  the  model  parameters  for  that  clay.  Some  data  close  to 
failure  have  been  ignored  in  selecting  the  design  curves  because  they 
are  not  consistent  with  the  rest  of  the  data.  This  inconsistency  may  be 
due  to  experimental  difficulties  in  capturing  failure  states  in  stress- 
controlled  testing  devices. 

The  experimental  tests  had  been  conducted  on  cylindrical  samples  in 
a  torsional  shear  testing  device.  All  samples  had  been  uni sotropical ly 
consolidated  to  insitu  conditions  and  all  the  tests  were  stress  con¬ 
trolled  and  performed  under  undrained  conditions. 

Figure  37  shows  model  predictions  for  a  torsional  shear  test  in 
which  the  major  principal  stress  was  inclined  at  0  =  15  degrees,  rela¬ 
tive  to  the  vertical  axis  of  the  soil  specimen  (Ref  40).  Figure  37  also 
siows  a  comparison  between  predicted  and  observed  behavior  of  the  soil 
in  these  tests.  Note  that  all  the  model  predictions  agree  well  with  the 
experimental  test  results  (Ref  40). 
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LABORATORY  CENTRIFUGE  TESTS 


Analysis  of  Model  Tests 


The  objective  of  this  work  was  directed  at  examining  the  ptedictive 
capabilities  of  the  Princeton  University  Effective  Stress  Soil  model. 
Of  particular  interest  are  the  validity  of  the  numerical  model  in 
capturing  the  generation  and  dissipation  of  excess  pore-water  pressures 
in  saturated  sand  deposits  during  and  after  earthquakes,  and  its  per¬ 
formance  in  dynamic  soil-structure  interaction  problems.  The  most 
appropriate  method  for  such  a  validation  study  would  be  to  field  data 
from  instrumented  prototype  situations.  However  such  a  study  is  pre¬ 
empted  by  the  nonava i 1 ab i 1 i ty  of  field  data.  In  the  absence  of  actual 
prototype  earthquake  field  data,  an  alternate  method  of  validation  is 
provided  by  analyzing  centrifugal  soil  model  tests  Although  imperfect 
in  many  respects,  it  is  felt  that  dynamic  centrifuge  soil  model  tests 
can  still  provide  a  valuable  data  base  for  calibration  of  numerical 
procedures  A  number  of  dynamic  centrifuge  soil  model  tests  have  been 
reported  in  the  literature  The  particular  tests  selected  for  this 
calibration  study  have  been  reported  in  References  ill,  ?? .  ?  1 .  and  ,”’4 
Monterey  and  the  L e i qhton- Buz zard  sands  were  used  in  these  tests  1  he 

basic  plan  developed  to  achieve  the  research  objective  is  summarized  as 
foil ows 
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A  summary  of  major  input  parameters  for  each  of  the  centrifuge  mode 
test  series  is  given  in  Table  2 

E  f  f  ec  1 1  ve- bt  ress  Soil  Model  Parameter's  The  ef fee t ’ ve- s t res s  s 
model  (Ret  Tl)  is  designed  such  that  all  its  required  pa’amete’s  .an  :• 
derived  from  the  results  of  standard  drained  compression  e*te"s  , 
t  r  i  a  *  ’  a  1  so’  ’  tests  data  ’he  ’'esu'-ts  from  the  i  onveH  •  ona  '  t  •  •  a  •  •  a 
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Table  2.  Summary  of  Analysis  Cases 
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Table  3.  Soil  Parameters 


Monterey  "0"  Sand  -  Or  =  40% 


's 

2.64  x 

103 

kg/m 

w 

1.00  x 

103 

kg/m 

w 

0.43 

k  = 

5  x  10‘ 

4  m/sec 

i‘  = 

34° 

c 

'e 

17.6° 

'1  = 

800.00 

1  = 

533.33 

1  = 

9.81  x 

104 

N/m2 

n  = 

0.50 

c 

1.20 

E  = 

0.60 

(mass  density,  solid  grains) 

(mass  density,  fluid  phase) 
(porosity) 

(permeabi 1 ity) 

(friction  angle;  compression  case) 
(friction  angle;  extension  case) 
(elastic  shear  modulus) 

(elastic  bulk  modulus) 

(reference  pressure) 

(power  exponent) 


Yie'd  Surface 
Number 


or 

m 

<WPl 

(H^j/p 

0.08937 

0.08937 

1181.00 

1581.00 

C. 18613 

0. 18610 

524.40 

I 

813.60 

0. 14217 

0.40800 

116.90 

210.50 

0  19995 

0.49770 

69. 15 

142.40 

26444 

0.63710 

23.46 

57.76 

;65  70 

0.81520 

3.  14 

10.08 

0  90910 

1.  11 

4.21 

■  '  ■  96990 

0.  52 

2.  17 

--  A  -%  w  **  /•  «  • 


Table  4  Soil  PsrjMtris 


Leigbton-Bui/ard  120/200  Sand  -  Or  =  55% 


ps 

= 

2  73  * 

103 

kg/m3 

(mass  density,  solid  grains) 

P» 

= 

1  00  * 

io3 

kg/m3 

(mass  density,  fluid  phase) 

w 

n 

= 

0  47 

-  3 

(porosity) 

k 

= 

2.5* 

10 

m/sec 

( permeab i H  ty ) 

= 

34  35° 

(friction  angle,  compression  cas 

= 

21. 14° 

(friction  angle,  extension  case) 

VP1 

= 

500. 

(elastic  shear  modulus) 

Bl/Pl 

= 

333. 33 

(elastic  bulk  modulus) 

P1 

= 

9.81  x 

io4 

N/m2 

(reference  pressure) 

n 

= 

0.  50 

(power  exponent) 

n 

— 

1.  30 

c 

nE 

= 

0.80 

Yield  Surface 
Number 


or 

m 

(HcVpi 

(H^j/p 

0. 03099 

0. 06582 

966.80 

1168. 00 

0.06261 

0. 14300 

357. 70 

459. 20 

0. 09681 

0.26590 

219.40 

337. 30 

0.08969 

0.32670 

162. 00 

265.40 

0. 15857 

0.50160 

42.03 

81.95 

0.23864 

0.63800 

17.47 

41. 31 

0. 27796 

0.76120 

7.60 

21. 99 

0. 32030 

0. 9023 

3.  78 

13.21 

0.32580 

0.9930 

1.  17 

4.47 

Figure  38.  Undeformed  mesh  soil  column  (Monterey  '0  sand). 

pore-fluid  is  allowed  to  take  place  through  the  rigid  bottom  boundary 

(nor  the  side  boundaries)  and  the  ground  shaking  is  applied  as  a 

horizontal  sinusoidal  input  acceleration  at  the  bottom  boundary  nodes. 

4 

The  permeability  (Table  3)  used  in  the  simulation  k  =  5  x  10  x 
100  m/sec  in  order  to  correct  for  the  fact  that  in  the  centrifuge  (at 
100  gs)  diffusion  of  pore-water  takes  place  100  times  faster  than  in  the 
corresponding  prototype. 

Figure  39  shows  the  computed  horizontal  acceleration  t  ^  histories 
at  the  bottom  (Figure  39(a))  and  at  the  top  (Figure  39(b))  of  the  soil 
column  for  an  input  horizontal  base  acceleration  with  an  amplitude 
0.285  gs  and  a  frequency  of  5  Hz.  The  results  for  25  cycles  of  loading 
(5  seconds  of  shaking)  are  reported  in  Figure  39.  Note  the  strong 
modification  of  the  signal  computed  at  the  surface  as  a  result  of  its 
passage  through  the  saturated  soil  deposit.  The  computed  maximum  sur¬ 
face  acceleration  compares  favorably  with  the  recorded  value  in  the 
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Figure  40  shows  the  computed  vertical  motion  at  the  surface.  As  a 
result  of  the  shaking,  excess  pore-water  pressures  build  up  and  partly 
dissipate  in  the  soil  column.  These  in  turn  generalize  vertical  motions 
(via  volumetric  strains  in  the  soil  skeleton).  Although  small,  the 
computed  vertical  acceleration  (Figure  40(b))  at  the  top  (0.0377  g)  is 
not  negligible  and  is  about  one-half  the  horizontal  peak  acceleration. 
Figure  40(a)  shows  the  resulting  settlement  of  the  soil  column  that 
accumulated  during  and  after  the  shaking  as  excess  pore-water  pressure 
is  being  dissipated  from  the  column.  The  computed  ultimate  settle¬ 
ment  =  24  cm. 

Figure  41  shows  the  computed  excess  pore-water  pressures,  vertical 
effective  stress,  and  shear  stress  time  histories  at  various  depths. 
The  plots  have  been  normalized  by  dividing  the  quantities  of  interest  by 
the  initial  vertical  effective  stress.  In  Figure  41(a),  h  =  0.51  meter 
and  in  Figure  41(b),  h  =  10.66  meters  where  h  =  distance  from  the  bottom 
boundary  (close  to  points  B  and  A  in  Reference  21.  Note  the  computed 
increases  in  pore-water  pressures.  In  Figure  41(a),  Au/Vq  =  0.69  and 
in  Figure  41(b),  Au/a^  =  0.85  compares  most  favorably  with  recorded 
values  ir.  the  test  (0.74  and  0.86,  respectively).  Note  the  diffusion 
taking  place  after  the  base  motion  has  stopped. 


Leighton-Buzzard  Soil  Column  Test 

The  test  procedures  and  test  results  are  reported  in  Reference  23. 
The  sand  was  rained  in  water,  in  a  stacked-ring  apparatus.  The  model 
was  tested  on  a  centrifuge  at  a  centrifugal  acceleration  of  35.5  gs,  and 
subjected  to  a  decaying  sinusoidal  base  acceleration.  The  corresponding 
prototype  situation  is  analyzed.  The  particular  test  selected  for 
analysis  is  referred  to  as  PL-3A  in  Reference  23. 

As  in  the  previous  section,  the  test  is  intended  to  simulate  free- 
field  conditions  in  a  horizontally  layered  soil  deposit,  and  the  same 
analysis  procedure  was  used.  Figure  42  shows  the  finite  element  mesh. 
Ten  equally  spaced  elements  are  used  to  simulate  the  10.8  meters  of  the 
sand  column  (12-inch  model  at  35.5  gs).  The  water  table  is  located  at 
the  ground  surface.  No  drainage  of  the  pore  fluid  is  allowed  to  take 
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Figure  42.  Undeformed  mesh  soil  column  (Leighton  Buzzard  sand). 

place  through  the  rigid  bottom  boundary  or  the  side  boundaries,  and  the 
ground  shaking  is  applied  as  a  horizontal  input  acceleration  at  the 
bottom  boundary  nodes.  The  permeability  (Table  4)  used  in  the  simula- 
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tion  is  k  =  2.5  x  10  x  35.5  m/sec  to  correct  for  the  fact  that  in  the 
centrifuge  (at  35.5  gs)  diffusion  of  the  pore  water  takes  place  35.5 
times  faster  than  in  the  corresponding  prototype. 

Figure  43  shows  the  computed  horizontal  acceleration  time  histories 
at  the  bottom  (Figure  42(a))  and  at  the  top  (Figure  42(b))  of  the  soil 
column  for  an  input  horizontal  base  acceleration  with  a  peak  amplitude 
of  0.20  g  and  a  frequency  of  1.72  Hz.  The  results  for  15  cycles  of 
loading  (9  seconds  of  shaking)  are  reported  in  Figure  42.  Note  that  as 
a  result  of  the  massive  liquefaction  occurring  in  the  underlying  strata, 
the  motion  is  not  transmitted  to  the  surface  after  about  2  seconds  of 
shaking  (1.16  cycles  of  loading)  as  observed  in  the  test.  Surface 
accelerations  were  not  measured  in  the  test,  which  make  direct  com¬ 
parison  with  the  predicted  time  history  shown  in  Figure  43(b) 
impossible. 

Figure  44  shows  the  computed  vertical  motion  at  the  surface.  Very 
little  iiffusion  takes  place  and  the  resulting  vertical  motions  and 
settlement  are  very  small  (7.3  cm). 


Figure  43.  Horizontal  accelerations  at  bottom  and  top  of  column 
(Lei ghton-Buzzard  sand). 
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gure  44.  Vertical  displacements  and  accelerations  at  top  of  col 
(Leighton-Buzzard  sand). 


Mgure  45  shows  the  computed  excess  pore-water  pressure,  vertical 
effective  stress,  and  shear  stress  time  histories  at  various  depths 
The  plots  have  been  normalised  by  dividing  the  quantities  of  interest  by 


the  initial  vertical  effective  stress  In  Figure  45(a).  h  -  H/'3  and  in 
Figure  45(b),  h  =  2 H/3  from  the  bottom  boundary,  where  H  -  10  8  meters  - 
height  of  soil  column  Note  the  rapid  increase  (in  almost  one  cycle  of 
loading)  in  pore-water  pressures  and  associated  liquefaction  of  the  soil 
column,  which  result  in  no  shear  being  transmitted  through  the  column 
The  predicted  massive  liquefaction  of  the  column  was  actually  recorded 
at  the  test  of  the  same  time,  and  resulted  in  failure  of  the  stack-ring 
apparatus . 

Brass  Footing  Test 

The  test  procedures  and  test  results  are  reported  in  Reference  24. 
The  soil  was  Lei ghton-Buzzard  120/200  sand.  The  soil  deposit  was  placed 
in  a  stacked-ring  apparatus  by  pluviating  the  sand  in  layers  into  water 
and  then  rodding  to  achieve  the  desired  density.  A  brass  footing  with  a 
diameter  of  113  mm  was  placed  on  top  of  the  saturated  sand  deposit 
(height  =  151  mm,  diameter  =  406  mm)  and  tested  on  a  centrifuge  at  a 
centrifugal  acceleration  of  80  gs.  The  deposit  was  then  subjected  to 
sinusoidal  base  acceleration.  The  corresponding  prototype  situation  was 
analyzed. 

Figure  46  shows  the  finite-element  mesh  used  for  analysis,  pore- 
pressures,  and  vertical  stresses  at  measured  points  and  a  comparison  of 
measured  and  computed  results.  The  soil  is  discretized  by  using 
240  elements  and  the  brass  footing  by  using  two  rows  of  10  elements 
each.  The  soil  parameters  are  given  in  Table  4,  and  in  the  analysis 
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k  =  2.5  x  10  x  80  m/sec  to  properly  scale  diffusion  time.  The  brass 

3  3 

footing  is  modeled  as  a  one-phase  elastic  solid  with  p  =  8.5  x  10  /kg/m 

12  2 

(mass  density);  E  =  10  N/m  (Young's  modulus),  and  v  =  0.0  (Poisson's 

ratio).  A  static  pressure  is  applied  to  the  top  of  the  footing  to 

5  2 

achieve  a  static  bearing  pressure  of  1.30  x  10  N/m  in  the  test  (at 


The  water  table  is  located  at  the  ground  surface.  Drainage  of  the 
pore  fluid  is  not  allowed  to  take  place  through  the  rigid  bottom 

boundary  or  the  lateral  side  boundaries.  A  ground  shaking  is  applied  as 
a  horizontal  sinusoidal  input  acceleration  at  the  bottom  boundary  nodes, 
with  a  maximum  acceleration  =  0.17  g  and  a  frequency  of  1  Hz  ( o' 

10  seconds  (10  cycles). 

The  stacked-ring  apparatus  is  used  to  simulate  free-field  cona¬ 
tions,  therefore,  the  same  procedure  as  used  previously  in  the  so- 
column  test  simulations  was  used.  Specifically,  the  same  eqoat -■ 

number  was  assigned  to  each  nodal  degree  of  freedom  on  the  s am* 
horizontal  plane  for  the  two  side  boundaries. 

Figure  47  shows  the  computed  horizontal  acceleration  t'me  -  - •  ■  .  . 
at  the  bottom  (Figure  47(a),  node  131)  of  the  soil  depos  - t  a-''  <■  •  . 

top  (Figure  47(b),  node  285)  of  the  brass  footing  on  the  e  * - 

Note  that  as  a  result  of  the  dynamic  soil-structure  • 

recorded  motion  on  the  footing  is  different  from  the  ■np,  *  -n  • 
base. 

Figure  48  shows  the  computed  vertical  acce'e-a'  •  •*. 
at  the  left  corner  (Figure  48(a).  node  275>  *> 

node  285),  and  right  corner  (Figure-  48(  c  '  . . 

Note  that  as  a  result  of  the  horizcmta'  have  •  .. 
are  imparted  to  the  footing  f  'qu^e  43 
vertical  displacements  at  the  e*r  •  * 

(Figure  49(b)),  and  right  rc-ie’  ;  -  •• 

test,  the  settlement  inc-easec 
the  ‘  :king  continued  l. 
after  the  shaking  st  -* 

HOWe  ,  er  ,  the  '  I'fTll  *  *■-  '  11  4  .  t.  rr 
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different  preparation  techniques,  resulted  in  differences  in  void  ratios 
for  the  respective  soils.  This  factor  effects  the  soils  ultimate 
strength  and  moduli.  The  computed  results  showed  stress  levels  on  the 
relatively  flat  yield  portion  of  the  stress-strain  diagram  where  small 
increases  in  loading  results  in  large  increases  in  displacement.  The 
actual  centrifuge  loading  conditions  remained  slightly  below  those  which 
would  generate  "flow-type"  strains  and  large  displacements.  In  addi¬ 
tion,  the  authors  (Ref  24)  acknowledged  difficulties  in  precisely 
controlling  the  properties  of  the  centrifuge  soil  deposit. 

These  factors  combined  to  generate  computer  simulation  which, 
though  only  moderately  weaker  in  ultimate  strength,  caused  much  larger 
strains,  hence,  displacements  to  be  generated  to  reach  an  equilibrium 
state  of  stress. 

Figure  50(a,  b,  and  c)  shows  the  deformed  mesh  at  t  =  5  seconds, 
10  seconds,  and  15  seconds. 

Figure  51  shows  the  contour  of  the  pore-water  pressures  at 
t  =  0  seconds,  10  seconds,  and  15  seconds.  Note  that  as  observed  in 
the  test,  in  the  "free-field"  close  to  the  sides,  the  pore-water  pres¬ 
sure  rises  quickly  (Figures  51(b), (c)).  Directly  under  the  structure, 
the  pore-water  pressure  increase  is  slower  and  always  remains  smaller 
than  the  pore-water  pressure  in  the  free-field  at  the  same  elevation 
(Figures  51(b),  (c)).  Immediately  following  the  shaking,  the  excess 
pore-water  pressures  dissipate  rapidly  and  reach  their  steady  state 
conditions  5  seconds  after  the  end  of  shaking  (Figure  51(d)).  This 
condition  was  further  illustrated  in  Figure  46,  which  shows  time 
histories  for  the  vertical  effective  stress  and  excess  pore-water 
pressures  for  the  points  shown  on  the  mesh  (Figure  45(a)). 

Retaining  Wall  Test 

The  test  procedures  and  test  results  are  reported  in  Reference  22. 
Loighton-Buzzard  sand  was  used  and  was  poured  dry  from  a  hopper  behind 
the  wall.  Density  was  adjusted  by  altering  the  rate  of  flow  and  height 
of  the  drop.  The  loose  backfill  case  (D^  =  55%)  was  selected  for  this 
analysis.  A  reinforced  micro-concrete  wall  model,  175  mm  high  with  a 


stem  thickness  of  15  mm  was  bolted  rigidly  to  the  test  container  with 
the  sand  backfill  placed  behind  it.  The  model  was  then  placed  on  board 
a  centrifuge  and  "spun  up"  to  a  centrifugal  acceleration  of  80  gs.  The 
container  was  then  subjected  to  a  sinusoidal  input  acceleration  motion 
perpendicular  to  the  plane  of  the  wall.  The  corresponding  prototype 
situation  is  analyzed  hereafter. 

Figure  52  shows  the  finite  element  mesh  used  for  the  analysis. 
The  backfill  is  discretized  by  using  280  elements.  The  soil  para¬ 
meters  are  given  in  Table  4.  The  retaining  wall  is  modelled  by  using 

3 

12  linear  beam  elements  with  p  =  2. 76  10  kg/m  (mass  per  unit  length). 
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Figure  52.  Undeformed  mesh  used  for  retaining  wall  test 
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E  =  1.70  x  10  N/m  (Young's  modulus),  G  =  6.54  x  103  N/in  (shear 

-14  2 

modulus),  I  =  1.44  10  m  (bending  moment  of  inertia);  and  A  =  1.20  m 

(cross-section  area).  The  soil  and  wall  subdomains  are  interfaced  by 

using  14  slide  line  elements,  with  K  =  10"1  (penalty  parameters).  Also, 

the  interface  between  the  wall  and  the  backfill  is  assumed  frictionless. 

The  bottom  boundary  is  assumed  rigid  and  the  side  boundary  smooth 
in  the  vertical  direction. 

The  computed  wall  crest  deflection  under  gravity  load  is  4.55  cm, 
which  compares  favorably  with  the  measured  test  value  (4.7  cm). 

In  order  to  simulate  the  dynamic  test  conditions,  the  same  sinu¬ 
soidal  horizontal  acceleration  was  imposed  on  the  bottom  boundary  nodes 
and  the  soil  right  side  boundary,  with  a  maximum  acceleration  =  0.20  g 
and  a  frequency  of  0.75  Hz  for  16  seconds  (12  cycles). 

Figure  53  shows  the  computed  horizontal  acceleration  time  histories 
at  the  bottom  (Figure  53(a))  and  at  the  top  (Figure  53(b))  of  the  wall. 
Note  that  as  the  result  of  the  interaction  with  the  nonlinear  soil  mass, 
the  computed  response  at  the  top  of  the  wall  is  amplified  and  exhibits 
superharmonics  (typical  of  nonlinear  systems).  Figure  54(a)  shows  the 
computed  horizontal  displacement  at  the  top  of  the  wall.  The  amplitude 
of  alternating  crest  deflection  is  3.725  cm  which  compares  most  favor¬ 
ably  with  the  recorded  amplitude  (0.50  x  80  =  4  cm)  in  the  test 
(Ref  22). 

Figure  54(b)  shows  the  computed  vertical  displacement  of  the  top 
soil  element  at  the  interface  with  the  wall.  Although  sinusoidal  in 
shape,  this  motion  occurs  with  a  much  longer  period  (about  12.5  seconds) 
than  the  input  horizontal  motion  (1.33  seconds). 

Figure  55  shows  the  computed  vertical,  horizontal,  and  shear  stress 
time  histories,  normalized  by  dividing  them  by  the  initial  vertical 
stress,  close  to  the  base  of  the  wall  (Figure  55(a),  element  A  in 
Figure  52)  and  at  midheight  (Figure  55(b),  element  B  in  Figure  52).  At 
midheight,  initially,  an  active  state  of  stress  prevails  and  K  =  0.40. 
As  a  result  of  the  stress  concentration  at  the  corner,  initially, 
K  =  0.675  at  the  bottom,  closer  to  a  passive  state.  As  a  result  of  the 
shaking,  the  lateral  stress  fluctuates  and  exhibits  a  net  decrease  at 
the  bottom  of  the  wall  (K  ~  0.35  after  12  cycles),  but  stays  about  the 
same  at  midheight. 


98 


Time  (set)  (x  10*) 
(  j) 


1  mu-  I  sve  m\  1  0  '  ) 
<hi 


Figure  53.  Input  horizontal  accelerations  and  resulting  accelerations 
at  top  of  retaining  wall. 
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MB 


Storage  Tank  Test 


The  centrifuge  model  consisted  of  a  storage  tank  placed  over  a  soft 
foundation  strata.  The  tank  model  was  constructed  of  rolled  aluminum 
plate,  which  formed  the  walls.  The  base  consisted  of  a  flexible,  rubber 
membrane  to  approximate  the  conditions  of  a  prototype  tank's  flexible 
base.  A  cross  section  of  the  model  package  is  presented  in  Figure  56. 
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Figure  56.  Model  tank  package 


The  model  was  then  mounted  on  a  centrifuge  and  accelerated  to 
60  gs.  The  model  was  maintained  at  this  acceleration  for  6,000  seconds 
to  allow  for  consolidation  to  take  place.  Once  hydrostatic  conditions 
were  reached,  the  simulation  of  the  tank  loading  and  unloading  began, 
Figure  56.  Pore  pressures  and  surface  displacements  were  recorded 
through  an  automated  data  acquisition  system.  More  detailed  information 
can  be  found  in  Reference  44. 

The  material  model  properties  for  the  kaolinite  simulation  were 
derived  from  undrained  triaxial  compression  and  extension  tests  pre¬ 
sented  above.  The  results  of  the  material  model  fitting  process  are 
shown  in  Figure  35.  The  triaxial  test  data  for  the  extension  test  was 
scaled  from  it's  original  confining  pressure  to  be  at  the  same  confining 
pressure  as  the  compression  test  as  required  for  input  to  Program  MUD. 
The  effect  of  this  scaling  on  the  final  solution  was  minimal. 

The  material  model  properties  for  the  layers  of  Monterey  "0"  sand 
were  constructed  from  test  data  presented  above.  Modifications  were 
made  to  the  critical  stress  ratio,  parameters  to  reflect  a  somewhat 
higher  relative  density,  which  was  back  calculated  from  the  weight  of 
the  sand  used  in  this  test.  Triaxial  test  data  were  not  available  for 
the  sand  used  in  this  test. 

The  finite  element  model  consisted  of  16  sand  and  48  clay  elements 
as  shown  in  Figure  57.  The  tank  loading  and  unloading  was  simulated  by 
a  uniform  pressure  distribution  over  the  inner  four  elements  as  shown  in 
Figure  57.  To  properly  initiate  the  stresses  and  consolidation  char¬ 
acteristics  of  the  model,  the  problem  was  started  at  an  initial  gravity 
of  zero,  and  then  increased  to  the  required  60  level  over  a  time  period 
of  250  seconds.  The  model  was  then  maintained  at  this  acceleration 
level  throughout  the  test.  After  the  appropriate  time  period  for 
consol idation  (6,000  seconds)  the  loading-unloading  simulation  began. 
The  loadingunloading  simulation  consisted  of  three  filling  and  two 
releasing  steps.  The  entire  trace  of  the  loading  is  shown  in  Figure  58. 
Each  loading  increment  consisted  of  a  125-second  load  step  and 
1,875  seconds  of  dissipation  to  allow  for  consolidation  time. 
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Figure  57.  Undeformed  mesh  used  for  model  storage  tank  test 
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Figure  58.  Pressure  loading  sequence. 


Because  this  is  a  layered  system,  which  is  similar  to  those  found 
in  nature,  it  is  important  to  correctly  capture  the  effects  of  all  of 
the  layers  on  the  system  response.  The  displacement  traces  shown  in 
Figures  59  and  60  show  the  instantaneous  settlement  of  the  sand.  Dis¬ 
placement  characterization  overall  is  in  good  agreement  while  the  actual 
measurements  are  acceptable.  There  was  some  difficulty  in  the  degree  of 
rebound  that  the  numerical  model  exhibited.  The  errors  here  are  caused 
by  the  difficulty  in  determining  the  correct  input  parameters  for  the 
material  model  as  well  as  soil  deposit  inconsistencies  in  placement  and 
uniformity  acknowledged  by  the  original  investigators.  The  pore- 
pressure  traces  (Figure  61)  are  in  excellent  agreement  with  those 
measured  during  the  test.  Dissipation  is  represented  with  excellent 
accuracy  in  all  cases  demonstrating  Program  DYNAFLOW's  capability  in 
modeling  threedimensional  consolidation  problems. 

Comparing  Program  DYNAFLOW's  prediction  and  those  of  the  finite 
element  code  used  by  the  original  investigators  (Reference  45)  show  the 
Program  DYNAFLOW  analysis  has  greater  capabilities.  This  is  primarily 
due  to  the  ability  of  the  Princeton  Effective  Stress  Soil  Model  to 
characterize  both  sands  and  clays,  as  opposed  to  modeling  the  sand  as  a 
quasi-elastic  soil  as  was  done  in  Reference  45.  The  prototype,  as 
simulated  in  this  configuration,  was  not  highly  sensitive  to  this 
assumption  of  elasticity,  as  the  problem  is  largely  dependent  on  only 
one  parameter--the  consolidation  coefficient.  However,  further  attempts 
at  modeling  more  complex  and  varied  systems  will  require  the  ability  to 
model  all  materials  present  in  the  system  as  is  done  in  the  Princeton 
Effective  Stress  Soil  Model. 

DYNAFLOW  ADINA  COMPARISON 

To  evaluate  the  relative  computational  efficiency  and  algorithm 
accuracy  of  Program  DYNAFLOW,  several  comparitive  total  stress  analyses 
were  run  using  DYNAFLOW  and  a  well  recognized  non-linear  finite  element 
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Figure  59.  Vertical  displacement  along  centerline. 
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Figure  60.  Free-field  displacements. 
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Figure  61.  Pore-pressure  traces. 
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DYNAF LOW- ADINA  comparison  load  versus  displacement,  case  I. 


Program,  ADINA.  Because  ADINA  has  only  total  stress  capability,  a 
constitutive  model  other  than  the  Princeton  Effective  Stress  Soil  Model 
had  to  be  selected.  The  Von  Mieses  model  is  included  in  both  computer 
programs  and  was  judged  best  for  a  straightforward  comparison. 

The  problem  selected  for  the  analyses  consisted  of  a  pressure 
loading  over  a  portion  of  a  semi-infinite  elastoplastic  half  space.  The 
problem  is  analagous  to  a  flexible  mat  foundation  on  which  the  distri¬ 
buted  load  is  increased  incrementally  until  bearing  capacity  failure  is 
reached. 
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The  theoretical  solution  for  limit  load  is  given  by  the  equation: 


P  =  (2  +  n)  C 

f 

where:  p^.  =  Failure  load  per  unit  area 

c  =  Cohesion,  or  maximum  shear  strength  of  the  material 
n  =  Pi 

The  analyses  were  run  in  three  configurations  to  determine  some  of 
the  sensitivities  to  changes  in  model  parameters.  The  relative  compu¬ 
tation  times,  material  parameters,  and  predicted  failure  loads  are: 


Analysis 

E 

C 

V 

(theo¬ 
retical  ) 

Program 

Pf 

(computed) 

cpu 

(seconds) 

I 

1000.0 
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II 

1000.0 
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III 

5000.0 

24.3 
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Adi  na 
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418 

E  =  Young' s  Modulus 
C  =  Cohesion 
v  =  Pc's  son's  Ratio 


The  resulting  load-deformation  curves  for  cases  I  and  II  are 
presented  in  Figures  62  and  63.  Note  in  Figure  62  an  additional 
analysis  using  Program  ADINA  is  included  with  a  plastic  modulus  of  1.0 
which  adds  a  significant  amount  to  the  total  stiffness  of  the  foundation 
system.  In  Figure  63  the  solutions  seem  to  change  their  respective 
"soft  versus  stiff"  relationship  present  in  Figure  62. 

This  case  is  a  representation  of  response  for  a  semi-incompressible 
material  (Poisson's  ratio  =  0.49).  The  stiffened  response  in  the 
program  ADINA  load-displacement  history  is  a  result  of  the  locking 
phenomona  generated  from  incompressible  elements.  As  Poisson's  ratio 
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Figure  63.  DYNAFLOW-ADINA  comparison  load  versus  displacement,  case  II. 

approaches  0.5  the  stiffening  becomes  more  pronounced.  The  solution 
using  Program  DYNAFLOW  involves  a  reduced  integration  on  the  volumetric 
stiffness  contribution,  and  prevents  the  locking  phenomona. 

Case  III  is  a  comparison  of  the  results  from  ADINA,  DYNAFLOW,  and  a 
solution  presented  in  Reference  44.  Poisson's  ratio  is  again  reduced  to 
0.3  and  the  differences  in  the  ADINA  and  DYNAFLOW  (Figure  64)  solutions 
become  more  pronounced.  The  differences  are  apparently  due  to  material 
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Figure  64.  DYNAFLOW-ADINA  comparison  load  versus  displacement, 
case  III. 

model  interpretation  in  the  respective  codes.  The  level  of  plasticity 
predicted  by  each  solution  is  the  same,  at  the  same  level  of  vertical 
deformation,  which  indicates  the  agreement  of  the  respective  strain 
level  computations  (Figure  64(b)).  The  stress  level  is  however, 
significantly  different  past  initial  yield  (Figure  64(a)). 


PILE  AND  PIER  FOUNDATIONS 


Load  Capacity 

This  discussion  considers  only  straight  piles  driven  into  homo¬ 
geneous  deposits  of  cohesionless  materials.  Experience  has  shown  that 
when  piles  are  driven  into  sand,  the  soil  near  the  pile  is  compacted  to 
a  distance  of  a  few  pile  diameters.  In  a  homogeneous  sand  tne  point 
resistance  and  average  skin  friction  increase  with  depth  of  penetration 
up  to  a  critical  depth.  Beyond  the  critical  depth,  the  point  resistance 
and  skin  friction  remain  almost  constant;  generally  this  is  caused  by 
soil  compressibility,  crushing,  and  arching.  The  empirical  approach  to 
prediction  of  pile  behavior  has  proven  more  satisfactory  than  an 
analytical  bearing  capacity  approach. 

The  load  transfer  mechanism  between  the  pile  and  the  surrounding 
soil  governs  the  behavior  of  the  pile.  The  design  of  a  pile  requires 
determination  of  the  proportion  of  load  transferred  to  the  soil  by 
adhesion  and  friction  between  the  pile  and  the  soil  and  that  transferred 
by  end  bearing.  This  is  influenced  by  the  flexibility  of  the  pile,  the 
stiffness  of  the  soil,  and  the  nature  of  the  transfer  mechanisms  between 
pile  surface  and  soil. 

Through  use  of  instrumented  field  and  model  studies,  the  design  of 
piles  and  piers  has  been  better  understood.  Figure  65  illustrates  use 
of  pile-driving  resistance  formulas  to  estimate  individual  allowable 
pile  loads.  Figure  66  illustrates  the  calculation  of  ultimate  load 
capacity  of  piles  in  cohesive  soils,  and  Figures  67  and  68  illustrate 
the  calculation  of  load  capacity  of  piles  in  cohesionless  soils.  How¬ 
ever,  these  conventional  ultimate  design  approaches  assume  the  simul¬ 
taneous  and  full  mobilization  of  pile  shear  resistance  and  base  bearing, 
which  is  not  well  founded.  Studies  have  shown  that  movement  of  a  pile 
must  be  present  to  mobilize  its  load-carrying  capacity. 

Reese  and  O'Neill  (Ref  47)  show  that  the  division  of  load  between 
sides  and  base  for  a  pier  on  stiff  clay  is  dependent  on  the  total  load¬ 
ing  (Figure  69).  The  first  two  load  increments  show  that  almost  all  of 
the  load  is  carried  by  friction;  as  the  load  increases,  more  is  carried 
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bearing  stratum  (see  example) 
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Figure  65.  Pile-driving  resistance  formulas  (Ref  45). 
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Figure  66.  Ultimate  load  capacity  of  piles  in  cohesive  soils  (Ref  45) 
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Figure  68.  Ultimate  load  capacity  of  driven  piles  in  noncohesive  soils 
(Ref  45). 
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Figure  69.  Load  distibution  for 
drilled  shaft  in  stiff  clay. 


by  end  bearing.  Above  80  tons,  any  additional  load  is  taken  by  end 
bearing.  Note  how  curve  C  simply  translates  to  form  curve  D,  indicating 
additional  load  is  taken  by  end  bearing.  With  the  information  in  Fig¬ 
ure  69,  a  typical  load  transfer  relationship  can  be  obtained  showing 
side  friction  (Figure  70).  From  the  slope  in  Figure  69,  it  is  evident 
that  frictional  effects  are  greatest  in  the  middle  (depth)  of  tne  pile., 
and  a  reduction  in  the  rate  of  load  transfer  occurs  in  the  lower  part, 
particularly  just  above  the  base.  The  distribution  of  frictional  forces 
depends  upon  the  soil  type.  Modifications  to  ultimate  strength  formula¬ 
tions  have  been  suggested  by  Reese  and  O’Neill  (Ref  47).  Factors 
include  an  effective  depth  concept  in  lateral  earth  pressure  calculation 
for  granular  soils  and  a  friction  reduction  factor  for  ultimate  resist¬ 
ance  in  cohesive  soils. 
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Figure  70.  Load  transfer  to  shaft 
segment  versus  movement. 


A  theoretical  solution  of  the  load  transfer  phenomenon  can  provide 
considerable  insight.  For  instance,  a  soil-structure  problem  would  be 
one  in  which  the  nature  and  mechanism  of  load  transfer  is  of  consid¬ 
erable  importance.  Solving  this  problem  must  focus  on  the  stress  strain 
characteristics  of  the  soil  and  the  behavior  of  the  interface  between 
pile  and  soil.  Using  a  linear-elastic  constitutive  model  would  be  a 
gross  simplification  of  the  real  material  properties.  The  Princeton 
University  Soil  Model  offers  an  opportunity  to  explore  nonlinear 
relationships. 

Several  approaches  have  been  used  to  model  the  contact  problem. 
Peterson  (Ref  48)  treats  the  two  contacting  surfaces  as  separate  and 
distinct  and  joins  them  mathematically  by  use  of  Lagrange  multipliers. 
However,  this  model  does  not  allow  preslip  deformation.  Herrman 
(Ref  49,  50,  and  51)  defines  three  behavior  modes:  nonslip,  slip,  and 
separation.  A  compatibility  model  combining  compatibility  and  equilib¬ 
rium  is  used.  The  compatibility  model  involves  linking  the  two  surfaces 
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with  fictitious  bond  springs,  and  the  frictional  forces  are  applied  as 
surface  tractions.  This  allows  a  preslip  deformation  that  would  not 
otherwise  be  computed. 

Another  approach  simply  links  two  nodes,  initially  overlapping  with 
nonlinear  springs.  The  procedure  is  simple  to  operate  but,  u. .fortu¬ 
nately,  may  produce  undesirable  numerical  characteristics  when  distinct 
rapid  changes  occur  in  stiffness. 

Prevost  has  formulated  a  contact  element  in  the  Program  DYNAFLOW. 
The  contact  element  may  be  used  to  impose  inequality  constraints  between 
nodes.  Either  friction  (i.e.,  "stick")  or  frictionless  (i.e.,  "slip") 
conditions  may  be  achieved. 

A  contact  element  is  defined  by  two  nodes:  a  spring  constant  or 
"penalty  parameter,"  k,  and  a  fixed  direction  vector,  +  d^  where  x^ 
is  the  initial  position  vector  and  d^  is  the  displacement  vector.  The 
contact  plane  passes  through  the  point  x^  is  the  displacement  vector. 
The  contact  plane  passes  through  the  point  x^  +  d^  and  is  perpendicular 
to  n  (Figure  71(a)).  The  contact/release  condition  is  defined  as 
fol lows: 

a  >  0  release 

a  §  0  contact 
where:  a  =  £  •  i n 

*  =  *B  *  sB  -  *A  -  sA 

The  quantity,  o,  is  a  measure  of  the  distance  between  Xg  +  dg  and 
the  contact  plane.  When  contact  is  noted,  a  contact  element  stiffness 
and  out-of-balance  force  are  added  to  the  global  equations. 

If  k  >  0  is  sufficiently  large,  the  point  Xg  +  dg  will  be  forced  to 
lie  (approximately)  on  the  contact  plane.  In  subsequent  steps,  only  the 
stiffness  is  assembled,  and  the  decision  to  remain  in  contact  or  to 
release  is  made  on  the  basis  of  the  sign  of  a,  as  above. 

For  interpreting  output,  the  contact  element  "displacement"  is 
defined  as  o,  and  the  "force"  is  given  by: 
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ka  if  o  <  0 

0  if  a  £  0 

In  addition  to  the  contact  element,  Prevost  (Reference  38)  has 
defined  a  slide  element.  The  slide-line  element  may  be  used  to  impose 
inequality  constraints  between  nodes.  Either  friction  or  frictionless 
conditions  may  be  achieved. 

A  slide- line  element  is  defined  by  three  nodes  and  a  spring  con¬ 
stant  or  penalty  parameter,  k.  The  connection  from  node  A  to  node  B 
defines  the  slide-line  direction,  and  node  C  is  the  contact  node 
(Figure  71(b)). 

The  projected  distance  of  node  C  to  node  A  onto  the  slide- line 
direction  is  denoted  by  a  and  is  given  by: 

a  =  •  %C  /  AB2  0  ^  a  1  1 

where:  •  =  the  dot  product  of  two  vectors.  The  direction  of  the  unit 

vector  n  to  the  slide-line  direction  is  given  by: 

n  =  (AS  x  A?)  x  AS  /  (A$  x  A6)  x  A§ 

where:  x  denotes  the  cross  product  of  two  vectors.  The  local  contact 

stiffness  matrix  k  is  given  by: 

(1  -  a)^  a(l  -  a)  -(1  -  a) 

k  =  k  a(l  -  a)  a  -or 

-(1  -  a)  -a  1 

where  the  rows  and  columns  are  arranged  such  that  the  first,  second,  and 
third  rows  (columns)  correspond  to  nodes  A,  B,  and  C,  respectively.  The 
contact/release  condition  is  defined  as  follows:  (1)  in  two  dimensions, 
if  0  s  a  i  and  Aff  •  n*  %  0,  contact  (otherwise,  release);  and  (2)  in 
three  dimensions,  if  0  S  a  i  1,  then  contact  (otherwise,  release). 

If  k  is  sufficiently  large,  the  point  C  will  be  forced  to  lie 
(approximately)  on  the  slide-line  AB.  In  subsequent  steps,  only  the 
contact  stiffness  is  assembled,  and  the  decision  to  remain  in  contact  or 
not  is  made  as  described  above. 
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Cyclic  Behavior  of  Piles 


Poulos  (Ref  52)  concludes  that  ultimate  load  capacity  and  cyclic 
stiffness  decrease  with  increasing  numbers  of  cycles  and  increasing 
cycle  load  level.  This  becomes  more  significant  when  the  cyclic  load 
approaches  one-half  the  static  ultimate  load.  The  cyclic  degradation 
appears  to  begin  at  the  top  of  the  pile  and  progresses  downward, 
resulting  in  a  gradual  transfer  of  load  to  the  lower  position  of  the 
pile.  The  crucial  factor  in  determining  the  amount  of  cyclic  degrada¬ 
tion  is  the  shear  strain  for  skin  friction. 

Poulos  (Ref  53)  conducted  a  limited  investigation  of  an  effective 
stress  approach  to  determine  pore  pressure  increases  with  cyclic  loading 
and  the  resulting  modulus  degradation  factors.  Figure  72  shows  a  compi¬ 
lation  of  observations  (Ref  52)  showing  degradation  as  a  function  of 
strain  ratio,  where: 

yss  =  (0.10  to  0.25)  ys 

and 

Ys  =  static  shear  strain  to  failure 

Poulos  (Ref  52)  reports  that  "one-way"  cyclic  loading  (i.e.,  cyclic 
loading  between  zero  minimum  load  and  a  specified  maximum)  produces 
"two-way"  loading  (i.e.,  loading  alternating  between  tension  and  com¬ 
pression  with  zero  as  a  mean  value).  He  points  out  that  degradation 
will  occur  at  different  rates  along  the  pile,  depending  on  local  stress 
level.  Even  in  initially  homogeneous  soil,  a  nonuniform  distribution  of 
soil  modulus  and  skin  friction  will  result  from  cyclic  loading  because 
nonuniformity  of  stress  dist,  ibution  occurs  along  the  pile.  Degradation 
occurs  in  the  ultimate  skin  friction  along  the  length  of  the  pile  and 
also  the  ultimate  base  resistance.  The  major  problem  in  a  cyclic 
response  is  determining  how  the  degradation  factors  vary  with  strain  and 
number  of  cycles. 

The  cyclic  shear  strain  in  the  soil  adjacent  to  the  pile,  Yr ,  can 
be  estimated  as: 
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cyclic  displacement  of  pile  at  a  point  on  the  pile  shaft 

pile  diameter 

In  [5  (1  -  ps)  L/d] 

1.0  homogeneous  infinitely  deep  soil  and  0.5  modulus 
increases  with  depth 

Poisson's  ratio 
embedder  pile  length 

Pi  .•  L  ood  (Wg) 

I0U  20C  300  400 

•  '  i  t  -  - 

‘i  4 

Horizontal  K  -  10"’ 

Vertical  K  =  10^ 

4- 

f>  t.  d. 

Figure  72.  Distribution  of  force  within  pile. 


The  cyclic  base  strain  can  be  estimated  from  the  above  assuming 
=  0.5,  L/d  =  100. 


where:  Pbc  =  cyclic  displacement  of  pile  base 
db  =  diameter  of  pile  base. 

It  is  important  to  note  that  the  above  is  based  on  the  assumptions 
of  elastic  behavior.  Typical  predictions  for  pile  settlement  take  the 
form  of  Figure  73.  Also  shown  in  Figure  73  is  the  ultimate  cyclic  load 
as  a  function  of  number  of  cycles.  The  settlement  problem  is  the  major 
concern  for  piles  in  cohesionless  materials. 

Friction  Between  Calcareous  Sand  and  Building  Materials 

Calcareous  sediments  have  proven  troublesome  to  offshore  facili¬ 
ties.  Piles  in  calcareous  sands  have  been  noted  to  penetrate  and  to  be 
extracted  with  much  less  effort  than  predicted  by  conventional  tech¬ 
niques.  A  research  program  was  previously  conducted  at  NCEL  (Ref  54)  in 
which  calcareous  sediments  were  collected  from  three  environments:  a 
deep-ocean  site,  a  shallow-ocean  site,  and  coral  line  sand  from  an  atoll 
beach.  The  coefficients  of  friction  of  these  sands  and  of  a  quartz  sand 
(used  as  a  standard)  were  measured  against  surfaces  of  rough  and  smooth 
steel  and  mortar.  Volume  changes  were  measured  as  a  function  of  sliding 
displacement. 

Experience  has  caused  engineers  to  reduce  pile  capacities  in  cal¬ 
careous  materials.  This  usually  results  in  load  capacity  reductions  to 
one-fourth  that  of  piles  in  normal  materials.  In  the  past,  it  was  not 
clearly  understood  why  calcareous  materials  exhibited  inferior  pile 
support.  Measured  angles  of  internal  friction  --  a  measure  of 
strength  -  are  high,  34  degrees  or  greater.  Part  of  the  problem  was 
identified  as  the  low  increase  in  soil  effective  stress  during  pile 
driving,  which  is  thought  to  result  from  a  crushing  or  collapse  of  a 
cemented  soil  structure  or  from  the  breakup  of  individual  carbonate 
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Figure  73.  Pile  degradation  with  cyclic 


grains.  This  lack  of  increase  in  effective  stress  results  in  a 
relatively  lower  shear  strength  in  the  soil  mass  surrounding  the  pile. 
Further,  it  was  thought  that  the  coefficients  of  friction  might  be 
substantially  less  between  calcareous  materials  and  piles.  Valent 
(Ref  54)  conducted  a  friction  test  program  using  a  modified  direct 
shear  test  machine  in  which  the  lower  half  of  the  apparatus  contained 
a  building  material  and  the  upper  half  a  soil  sample.  Table  5  sum¬ 
marizes  his  results. 

In  general,  the  results  show  that  the  low  friction  forces  in  cal¬ 
careous  sediments  are  not  the  result  of  low  achievable  friction  angles, 
because  calcareous  sands  are  comparable  to  other  sands;  low  friction 
development  then,  must  be  attributable  to  low  normal  force.  Deep-ocean 
sand  (foraminiferal  sand-silt)  exhibits  one  possible  cause  for  low 
developed  normal  force.  The  volume  change  during  testing  indicated  a 
considerable  volume  decrease  during  development  of  resisting  friction 
force,  probably  due  to  crushing  of  the  skeletal  structures  and  shell 
fragments.  Penetration  of  a  pile  in  such  a  material  would  crush  the 
hollow  shell  material  with  only  a  minimal  increase  in  effective  stress 
of  the  surrounding  material. 

Quoting  from  Valent  (Ref  54): 

"1.  The  calcareous  sediments  tested,  and  presumably 
calcareous  sediments  in  general,  develop  coefficients 
of  friction  against  steel  and  concrete  building 
materials  that  are  comparable  to  those  developed 
by  quartz-type  sands.  Thus,  the  possibility  of 
low  coefficients  of  friction  being  responsible  for 
the  observed  low  friction  forces  on  driven  piling 
and  other  penetrators  in  calcareous  materials  is 
ruled  out. 

2.  The  observed  large  volume  decreases  during 
shear  of  the  foraminiferal  sand-silt  are  probably 
responsible  for  the  low  developed  friction  forces 
in  these  hollow-shelled  materials.  Such  large 
volume  decreases  at  nonincreasing  normal  load  imply 
dens i fi cation  in  the  field  without  accompanying 
increases  in  normal  stress  on  the  penetrator  surface. 

3.  Low  developed  friction  forces  in  other  calcareous 
materials  may  arise  from  a  similar  mechanism  involving 
a  hypothesized  loose,  but  cemented,  structure  for 


Table  5.  Summary  of  Friction  Test  Results  (Ref  53) 


Test  No. 

Soil  Material 

a 

^peak 

a 

^residual 

Base  Material*3 

Sandc 

1 

Quartz  sand 

0. 67d 

0.54 

2 

Coralline  sand 

0.66 

0.56 

3 

Coralline  sand 

0.68 

0.  57 

4 

Oolitic  sand 

0.  77d 

0.61 

5 

Oolitic  sand 

0.81. 

0.62 

6 

Foram  sand- si  It 

0.  64d 

0.58 

Base  Material*3:  Smooth  Steel 

7 

Quartz  sand 

0.  27d 

0.19 

8 

Coralline  sand 

0.20 

0.17 

9 

Coralline  sand 

0.20 

0. 18e 

10 

Coralline  sand 

0.  21 

0.17 

11 

Oolitic  sand 

0.15d 

0. 13 

12 

Oolitic  sand 

0.32 

0.31 

13 

Foram  sand-silt 

0.40 

0.37 

Base  Material*3:  Rough  Steel 

14 

Quartz  sand 

0.60 

0.54 

15 

Coralline  sand 

0.63 

0.55 

16 

Oolitic  sand 

0.54 

0.51 

17 

Oolitic  sand 

0.  58, 

0.50 

18 

Foram  sand-silt 

0.66 

Base  Material*3:  Smooth  Concrete 

19 

Quartz  sand 

0.60 

0.54 

20 

Coralline  sand 

0.63 

0.56 

21 

Oolitic  sand 

0.59 

0.52 

22 

Oolitic  sand 

0.58, 

0.  54 

23 

Foram  sand-silt 

T 

0.67 

128 


Table  5.  Continued 


Test  No. 

Soil  Material 

a 

^peak 

a 

^residual 

Base  Material*1:  Rough  Concrete 

24 

Quartz  sand 

0.69 

m 

25 

Coralline  sand 

0.66 

26 

Oolitic  sand 

0.74 

■saH 

For  direct  shear  tests,  p  =  tan  0  where  0  =  angle  of  internal 
friction;  for  friction  tests,  p  =  tan  6  where  6  =  angle  of 
sliding  friction. 

Soil  in  bottom  shear  ring  for  direct  shear  tests,  or  building 
material  in  friction  tests. 

Base  material  same  as  soil  material  for  direct  shear  tests. 

! 

These  tests  run  with  mechanical  measurement  system; 
i.e.,  proving  ring  and  manual  recording  of  data. 

Low  value  for  p  reached  shortly  after  ^peak,  thereafter  p 
increased  with  displacement  to  end  of  test. 


No  peak  p  reached,  p  increasing  through  end  of  test. 


the  soil  material.  The  application  of  shear  stresses 
during  penetration  would  cause  collapse  of  this 
structure  to  a  denser,  but  still  loose,  arrangement. 11 


Piles  in  Calcareous  Sands 


Calcareous  sands,  as  discussed,  are  noted  for  loose  arrangement  of 
particles  lightly  cemented  to  form  a  structure  to  support  other  layers 
without  compacting.  However,  upon  shearing,  the  structure  is  destroyed, 
breaking  the  cement  bonds.  The  loose-grained  structure  then  compacts 
and  densities.  The  sand  (silica  sand)  is  a  very  loose  sand  (30%  rela¬ 
tive  density)  and,  as  such,  is  a  "manufactured"  sand  that  would  not 
occur  in  nature.  The  properties  of  this  material  might  be  expected  to 
represent  the  constitutive  behavior  of  a  calcareous  sand. 

To  analyze  a  pile's  load  capacity  using  the  effective  stress  soil 
mocel  ,  the  mesh  in  Figure  74  was  used.  This  mesh  simulates  a  pile 
already  in  place  (i.e.,  not  the  driving  of  the  pile).  Horizontal  and 
vertical  springs  were  used  to  join  the  pile  to  the  soil  field  to  simu¬ 
late  the  interface  and  allow  for  pile  movement. 

The  soil  properties  used  were  those  of  silica  sand.  The  material 
properties  simulate  the  condition  of  the  soil  after  placement  of  the 
pile.  No  desification  occurs,  since  the  pile  is  in  place  at  the  start 
of  the  analysis.  A  prediction  of  pile  capacity  was  made  by  using  con¬ 
ventional  techniques  found  in  Figures  67  and  68  adjusted  to  the  problem 
conditions  but  without  densi f ication.  Results  indicated  a  pile  resis¬ 
tance  in  normal  average  cohesionless  sand  of  at  least  1,500  kg  with 
about  150  kg  of  that  in  side  friction  and  the  remainder  in  end  bearing. 

Load  Settlement  of  the  Pile.  Figure  75(a)  shows  the  load  settle¬ 
ment  of  the  pile.  Several  values  were  tried  for  the  spring  constants. 
Results  show  the  pile  experienced  large  settlements  between  300  and  400 
kg.  Figure  72  shows  the  distribution  of  force  within  the  pile  for  three 
load  levels.  Note  that  most  of  the  increase  in  load  is  due  to  skin 
friction.  Figure  75(b)  shows  the  distribution  of  force  along  the  pile 
at  a  point  in  the  loading  having  a  force  of  300  kg  applied  at  the  top  of 
the  pile.  The  results  in  Figure  75(b)  show  about  one-half  to  two-thirds 
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of  the  load  is  taken  in  skin  friction.  The  skin  friction  portion  of  the 
load  agrees  with  the  predicted  values.,  however,  the  end  bearing  is 
substantially  less.  It  is  important  to  note  that  the  silica  sand  used 
had  a  relative  density  of  only  30%  which  is  so  loose  it  does  not  occur 
in  nature  but  is  "manufactured"  in  the  laboratory.  The  level  of  loading 
is  about  one-fourth  that  of  a  pile  in  normal  sand,  which  is  the  level 
expected  from  experience  in  calcareous  sands.  Figure  76  shows  contours 


of  principal  stress  around  the  pile  at  a  typical  load  step.  These  show 
typical  patterns,  as  expected;  the  level  of  stress  below  the  pile  tip 
would  be  on  the  last  surface  at  yield. 

Results  show  that  the  order  of  magnitude  of  the  pile  capacity  in 
calcareous  materials  is  predicted  correctly.  Further,  the  friction 
developed  on  the  sides  of  the  pile  is  at  the  expected  level,  not  sub¬ 
stantially  reduced  from  normal  sands.  The  amount  of  end  bearing  is 
slightly  dependent  on  the  vertical  spring  constant.  Use  of  nonlinear 
springs  would  have  allowed  slip  to  occur,  transferring  more  load  into 
end  bearing  when  exceeding  some  local  slip  level. 

Pile  Behavior  Under  Cyclic  Loading.  Next,  pile  behavior  under 
cyclic  loading  was  examined.  The  applied  loading  consisted  of  cyclic 
variation  of  the  vertical  load  on  the  pile.  In  one  case,  the  loading 
was  cycled  in  "one-way"  loading  (zero  minimum  and  specified  maximum)  and 
in  "two-way"  loading  (alternating  compression  and  tension).  The  first 
example  was  at  a  relatively  low  level  of  loading  for  the  model  pile 
(about  one-seventh  of  yield).  Figure  77  shows  the  stress  beneath  the 
pile  tip  for  one-way  loading  (0  to  50  kg),  Figure  78  shows  the  shear 
stress  contours  around  the  pile  (0  to  50  kg).  Figure  79  shows  the 
stress  beneath  the  pile  tip  for  two-way  loading  (50  to  50  kg),  the  shear 
stress  plot  is  similar  to  Figure  80.  Note  that  the  tip  forces  reduce 
with  loading,  transferring  more  to  friction  to  maintain  load  levels. 
Figure  80  shows  the  typical  contact  element  force,  a  measure  of  friction 
between  the  soil  and  pile.  Note  the  gradual  buildup  in  the  two-way 
loading  case. 

The  loading  was  increased  to  about  40%  of  yield.  Figure  81  shows 
the  distribution  of  force  within  the  pile  for  both  one-way  and  two-way 
loading  for  the  first  and  fifth  cycles.  Note  the  increase  in  friction 
loading  in  the  pile  with  two-way  loading,  showing  the  degradation  with 
cyclic  loading.  Note  also  that  for  both  cases  the  friction  increases 
more  in  the  lower  half  of  the  pile.  Figure  82  represents  the  typical 
contact  element  force,  again  showing  the  increase  in  friction  with  each 
cycle.  Figures  83  and  84  show  the  stress  beneath  the  pile  tip.  Both 
drop  off  with  each  cycle,  however,  the  two-way  loading  does  so  at  a 


faster  rate.  Figure  85  shows  the  pile  settlement  for  both  cases.  Note 
that  the  change  in  settlement  between  the  first  and  fifth  cycle  is  over 
four  times  greater  for  the  two-way  loading.  Figure  86  shows  the 
deformed  mesh  for  the  one-way  loading.,  the  mesh  is  similar  for  the 
two-way  loading. 

The  loading  was  increased  to  80%  of  pile  capacity,  and  the  one-way 
loading  case  was  tried.  This  load  exceeded  the  friction  capacity  so 
two-way  loading  could  not  be  used.  Figure  87  shows  the  stress  beneath 
the  pile  tip,  and  Figure  88  shows  the  deformed  mesh.  Failure  occurs 
between  load  steps  70  and  80.  Figure  89  shows  the  shear  stress  contours 
around  the  pile  at  load  steps  10  and  70.  More  load  capacity  is  trans¬ 
ferred  from  end  bearing  than  from  side  friction. 

The  cyclic  results  clearly  show  the  degradation  of  the  pile  under 
cyclic  loading.  Results  show  the  two-way  loading  causes  more  degrada¬ 
tion,  as  expected.  The  model  performed  well  in  predicting  the  qualita¬ 
tive  pile  behavior,  particularly  the  settlements  (Figures  75  and  85). 

The  pile  problems  were  repeated  for  the  undrained  case  in  which  the 
water  table  was  at  the  surface.  Figure  90(a)  shows  the  load  settlement 
curve.  The  original  load  increment  used  in  the  drained  problem  above 
was  too  large  for  the  undrained  problem  and  produced  an  instability  at 
yielding.  The  load  step  was  reduced,  and  the  solution  proceeded  in  a 
satisfactory  manner.  The  pile  capacity  was  about  250  kg,  lower  than  the 
350  kg  in  the  drained  case.  The  loading  was  cycled  in  one-way  and 
two-way  conditions.  Figure  90(b)  presents  the  distribution  of  force 
within  the  pile,  showing  the  increase  in  friction  with  cycling,  as  with 
the  drained  case.  Friction  is  slightly  greater  in  this  case. 
Figures  91  and  92  show  the  soil  stress  beneath  the  pile.  The  end 
bearing  stress  drops  off  with  cycling,  showing  the  shift  to  friction. 
The  two-way  loading  has  a  faster  degradation,  as  clearly  shown  in 
Figure  93,  which  shows  the  settlement.  The  two-way  loading  produces 
more  permanent  differential  settlement  after  the  fifth  cycle.  Figure  94 
shows  the  pore  pressure  in  the  soil  beneath  the  pile  tip.  The  pore 
pressure  rises  to  about  three-quarters  of  the  confining  stress  and  then 
drops  off  as  the  load  is  shifted  to  friction  along  the  pile.  Figure  95 
shows  the  friction  force.  Note  the  degradation  effects  with  cycling  as 
the  load  transfer  drops  off  with  each  cycle. 
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Figure  79.  Horizontal  and  vertical  stress  of  pile  in  soil 
field  -  -50  to  50  kg. 
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Figure  80.  Friction  force  in  pile  in  soil  field  -  0  to  50  kg 
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Figure  83.  Horizontal  and  vertical  stress  of  pile  in  soil  field 
0  to  150  kg. 
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Figure  85.  Horizontal  and  vertical  settlement  stress  of  pile  in  soil 
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Figure  86.  Deformed  mesh,  cyclic  load  -  0  to  50  kg. 
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Figure  87.  Horizontal  and  vertical  stress  of  pile  in  soil  field 
0  to  300  kg. 
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Figure  88.  Deformed  mesh,  cyclic  load  -  0  to  300  kg. 
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Figure  93.  Pile  settlement  under  cyclic  load  (with  water  table) 
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Figure  95.  Friction  force  of  pile  in  soil  field  (with  water  table). 


CONCLUSIONS/RECOMMENDATIONS 

An  efficient  finite  element  procedure  to  analyze  dynamic  transient 
phenomena  in  dry  and/or  fluid-saturated  porous  soil  media  is  presented. 
The  saturated  porous  medium  is  modeled  as  a  two-phase  system  consisting 
of  a  solid  and  a  fluid  phase.  Time  integration  of  the  resulting  semi¬ 
discrete  finite  element  equations  is  performed  by  using  an  implicit- 
explicit  algorithm.  In  order  to  remove  the  time  step  size  restriction 
associated  with  the  presence  of  the  stiff  fluid  in  the  mixture,  the 
fluid  contribution  to  the  equations  of  motion  is  always  treated 
implicitly.  The  procedure  allows  an  optimal  selection  of  the  time  step 


size  independently  of  the  fluid.  Depending  upon  the  particular  appli¬ 
cations  (seismic,  blast  loading,  etc.)  the  fluid  may  be  assumed 
incompressible  or  compressible. 

Accuracy  and  versatility  of  the  proposed  procedure  are  demonstrated 
by  applying  it  to  a  number  of  dynamic  soil  and  soil-structure  inter¬ 
actions  tests  performed  in  centrifuges.  The  hysteretic  stress-strain 
behavior  of  the  soil  skeleton  is  modeled  by  using  the  ef fecti ve-stress 
elastic-plastic  model  reported  in  Reference  31.  The  model  is  an 
extension  of  the  simple  multi-surface  ^‘Plasticity  theory  and  uses 
conical  yield  surfaces.  The  model  accurately  describes  observed  shear 
nonlinear  hysteretic  soil  behavior,  shear  stress- i nduced  anisotropy 
effects,  and  reflects  accurately  the  strong  dependency  of  the  shear 
dilatancy  on  the  effective  stress  ratio  in  granular  cohesionless  soils. 
The  model  is  applicable  to  general  multi-dimensional  stress-strain 
conditions,  and  its  parameters  can  be  derived  entirely  from  the  results 
of  conventional  triaxial  soil  tests. 

The  proposed  procedure  allows  a  complete  and  rational  analysis  of 
complex  dynamic  soil-structure  interaction  problems  including  important 
effects  such  as  hysteretic  nonlinear  effective  stress-strain  behavior, 
coupling  between  shear  and  volumetric  deformations,  coupled  effects  in 
saturated  two-phase  soil  systems. 

The  soil  model  has  been  validated  by  comparison  with  laboratory 
test  data,  predicting  drained  tests,  undrained  tests,  and  tests  at 
differing  stress  paths.  The  model  has  been  further  validated  by 
comparison  of  computed  results  to  results  obtained  in  centrifuge  model 
tests. 


REFERENCES 

1.  Biot,  M.A.,  "Theory  of  propagation  of  elastic  waves  in  a  fluid- 
saturated  porous  solid,"  Journal  of  Acoustical  Society  of  America, 
vo 1  28,  1956,  pp  168-191. 


155 


2.  Van  der  Kogel,  H.  ,  "Wave  propagation  in  staurated  porous  media," 
PhD  Thesis,  California  Institute  of  Technology,  Pasadena,  Calif.,  1977. 

3.  Nikolaeviskii ,  V.N.;  Vasniev,  K.S.;  Gorbunov,  A.T. ;  and  Zotov,  G.A. 
Mechanics  of  saturated  porous  media,  Nedra,  Moscow,  1970,  335  pages. 

4.  Garg,  S.K. ;  Brownell,  D.H.;  Pritchett,  J.W.  ;  and  Herrmann,  R.G., 
"Shock-wave  propagation  in  fluid-saturated  porous  media,"  Journal  of 
Applied  Physics,  vol  46,  no.  2,  1975,  pp  702713. 

5.  Garg,  S.K. ;  Nayfeh,  A.H.;  and  Good,  A.J.,  "Compressional  waves  in 
fluid-saturated  elastic  porous  media,"  Journal  of  Applied  Physics, 
vol  45,  no.  5,  1974,  pp  1968-1974. 

6.  Ghaboussi,  J.  and  Dikmen,  S.U. ,  "Liquefaction  analysis  of  hori¬ 
zontally  layered  sands,"  Journal  of  Geotechnical  Engineering  Division, 
American  Society  of  Civil  Engineers,  vol  104,  no.  FT3,  1978,  pp  341-356. 

7.  Liou,  C.P.;  Specter,  V.L.  ;  and  Richard,  F.E.,  "Numerical  model  for 
liquefaction,"  Journal  of  Geotechnical  Engineering  Division,  American 
Society  of  Civil  Engineers,  vol  103,  no.  GT6,  1977,  pp  589-606. 

8.  Bowen,  R.M.  and  Lockett,  R.R. ,  "Inertial  effects  in  poroelasticity," 
Journal  of  Applied  Mechanics,  American  Society  of  Mechanical  Engineering 
Transactions,  vol  50,  1983,  pp  334-342. 

9.  Mei ,  C.C.  and  Foda,  M.A. ,  "Wave- i nduced  responses  in  a  fluid-filled 
poro-elastic  solid  with  a  free  surface:  A  boundary  layer  theory,"  Geo¬ 
physics  Journal,  Royal  Astrological  Society,  vol  66,  1981,  pp  597-631. 

10.  Ghaboussi,  J.  and  Wilson,  E. L. ,  "Variational  formulation  of 
dynamics  of  fluid-saturated  porous  elastic  solids,"  Journal  of 
Engineering  Mechanics  Division,  American  Society  of  Civil  Engineers, 


11.  Lysmer,  J.  ,  "Analytical  procedures  in  soil  dynamics,"  in  Proceed¬ 
ings  of  the  American  Society  of  Civil  Engineers,  Geotechnical 
Engineering  Division,  Specialty  Conference  on  Earthquake  Engineering  and 
Soil  Dynamics,  Pasadena,  Calif.,  vol  3,  June  1978,  pp  1267-1316. 

12.  Prevost,  J.H. ,  "Effective  stress  analysis  of  seismic  site 
response,"  International  Journal  for  Numerical  Analytical  Methods  in 
Geomechanics,  Dec  1985. 

13.  Prevost,  J.H.,  "Wave  propagation  in  f luid-saturated  porous  media: 
An  efficient  finite  element  procedure,"  International  Journal  of  Soil 
Dynamics  in  Earthquake  Engineering,  Nov  1985. 

14.  Prevost,  J.H. ,  "Nonlinear  transient  phenomena  in  soil  media,"  in 
Mechanics  of  Engineering  Materials,  Chapter  26,  Editors:  C.S.  Desai  and 
R.H.  Gallagher.  New  York,  N.Y. ,  Wiley  and  Sons,  1984,  pp  515-533. 

15.  Prevost,  J.H.,  "Implicit-explicit  schemes  for  nonlinear  consolida¬ 
tion,"  Computer  Methods  in  Applied  Mechanics  Engineering,  vol  39,  no.  2, 
1983,  pp  225-239. 

16.  Prevost,  J.H. ,  "Nonlinear  transient  phenomena  in  saturated  porous 
media,"  Computer  Methods  in  Applied  Mechanics  Engineering,  vol  30,  1982, 
pp  3-18. 

17.  Zienkiewicz,  O.C.  and  Shionii,  T.  ,  "Dynamic  behavior  of  saturated 
porous  media;  the  generalized  biot  formulation  and  its  numerical 
solution,"  International  Journal  for  Numerical  Analytical  Methods  in 
Geomechanics,  vol  8,  no.  1,  1984,  pp  71-96. 

18.  Hughes,  T.J.R.  and  Liu,  W. K. ,  "Implicit-explicit  finite  elements  in 
transient  analysis:  Stability  theory,"  Journal  of  Applied  Mechanics, 
American  Society  of  Mechanical  Engineers  Transactions,  vol  45,  1978, 


19.  Hughes,  T.J.R.  and  Liu,  W.K.,  "Implicit-explicit  finite  elements  in 
transient  analysis:  Implementation  and  numerical  examples,"  Journal  of 
Applied  Mechanics,  American  Society  of  Mechanical  Engineers  Trans¬ 
actions,  vol  45,  1978,  pp  395-398. 

20.  Prevost,  J.H.,  "Study  of  the  capabilities  of  a  soil  model  to  predict 
axial  drained/undrained  responses  of  two  sands,"  Report  to  NCEL, 
June  1985. 

21.  Arulanandan,  K.  ;  Anandarajah,  A.;  and  Abghari ,  A.,  "Centrifugal 
modeling  of  soil  liquefaction  susceptibility,"  Journal  Geotechnical 
Engineering,  American  Society  of  Civil  Engineers,  vol  109,  no.  3, 
Mar  1983,  pp  281-300. 

22.  Bolton,  M.D.  and  Steedman,  R.S.,  "Centrifugal  testing  of  micro¬ 

concrete  retaining  walls  subjected  to  base  shaking,"  in  Proceedings  of 
the  Soil  Dynamics  and  Earthquakes  Engineering  Conference,  Southampton, 
U.K. ,  1982,  pp  311-329. 

23.  Lambe,  P.C.,  "Dynamic  centrifugal  modelling  of  a  horizontal  sand 

stratum,"  Sc.  D.  Thesis,  MIT,  Cambridge,  Mass.,  1981. 

24.  Whitman,  R.V.  and  Lambe,  P.C.,  "Liquefaction:  Consequences  for  a 

structure,"  in  Proceedings  of  the  Soil  Dynamics  and  Earthquakes 
Engineering  Conference,  Southampton,  U.K.,  1982,  pp  941-949. 

25.  Terzaghi ,  K.  ,  Theoretical  soil  mechanics.  New  York,  N.Y.,  Wiley 

and  Sons,  1943. 

26.  Bowen,  R.M.  ,  "Theory  of  mixtures,"  in  Continuum  Physics,  vol  III, 

edited  by  A.C.  Eringen,  Academic  Press,  1976,  pp  1-127. 


158 


27.  Eringen,  A.C.,  and  Ingram,  J.D.,  "A  continuum  theory  of  chemically 
reacting  media  I  and  II,"  International  Journal  of  Engineering  Science, 
vo 1  3,  1955,  pp  197-212,  and  vol  5,  1967,  pp  289-322. 


m 


28.  Green,  A.C.  and  Naghdi ,  P.M. ,  "A  dynamical  theory  of  interacting 
continua,"  International  Journal  of  Engineering  Science,  vol  3,  1965, 
pp  231-241. 

29.  Bowen,  R.M. ,  "Compressible  porous  media  models  by  use  of  the  theory 
of  mixtures,"  International  Journal  of  Engineering  Science,  vol  20, 
no.  6,  1982,  pp  697-735. 

30.  Hill,  R. ,  "A  general  theory  of  uniqueness  and  stability  in  elastic- 
plastic  solids,"  Journal  of  Mechanics  and  Physics  of  Solids,  vol  6, 
1958,  pp  236-249. 

31.  Provost,  J.H..  "A  simple  plasticity  theory  for  frictional  cohesion¬ 
less  soils,"  International  Journal  of  Soil  Dynamics  in  Earthquake 
Engineering,  vol  4,  no.  1,  1985,  pp  9-17. 

32.  Prevost,  J.H..  "Plasticity  theory  for  soil  stress-strain  behavior," 
Journal  of  Engineering  Mechanics  Division,  American  Society  of  Civil 
Engineers,  vol  104,  no.  EM5,  1978,  pp  1177-1194. 

33.  Zienkiewicz,  O.C.,  "The  finite  element  method,"  3rd  Edition, 
London,  McGraw-Hill,  1977. 

34.  Hughes,  T.J.R.  and  Belytschko,  T.  ,  "A  precis  of  developments  in 
computational  methods  for  transient  analysis,"  Journal  of  Applied 
Mechanics,  American  Society  of  Mechanical  Engineers  Transactions, 
vol  50,  no.  4b,  1983,  pp  1033-1041. 

35.  Hughes,  T.J.R. ;  Pister,  K.S.;  and  Taylor,  R.L.,  "Implicit-explicit 
finite  elements  in  nonlinear  transient  analysis,"  Computer  Methods  in 
Applied  Mechanics  Engineering,  vol  17/18,  1979,  pp  159-182. 


W  "a  *r  .*  -  •»  A  .%  J •  A  L*.  A  .  r  L  %  _  .  '  %  A  V  >  A  'V  A"  ■.  A  J 


36.  Newmark,  N.M. ,  "A  method  of  computation  for  structural  dynamics, 
Journal  of  Engineering  Mechanics  Division,  American  Society  of  Civil 
Engineers,  vol  85,  no.  EM3,  1959,  pp  67-94. 

37.  Hilber,  H.M. ,  "Analysis  and  design  of  numerical  integration  methods 
in  structural  dynamics,"  Report  No.  EERC  76-29,  Earthquake  Engineering 
Research  Center,  University  of  California,  Berkeley,  Calif.,  Nov  1976. 

38.  Prevost,  J.H..  "DYNAFLOW:  A  nonlinear  transient  finite  element 

analysis  program,"  Department  of  Civil  Engineering,  Princeton 

University,  1981.  (Last  revision  Aug  1985) 

39.  Atkin,  R.J.,  "Completeness  theorems  for  linearized  theories  of 
interacting  conti nua,"  Quarterly  Journal  of  Mechanics  and  Applied 
Mathematics,  vol  21,  1968,  pp  171-193. 

40.  Bowen,  R.M.  and  Reinicke,  K.M. ,  "Plane  progressive  waves  in  a 

binary  mixture  of  linear  elastic  materials,"  Journal  of  Applied 
Mechanics,  American  Society  of  Mechanical  Engineers  Transactions, 

vol  45,  1978,  pp  493-499. 

41.  Bowen,  R.M.  and  Chen,  P.J.,  "Waves  in  a  binary  mixture  of  linear 

elastic  materials,"  Journal  of  Medanique,  vol  14,  no.  2,  1975, 

pp  237-266. 

42.  Bathe,  K.J.  and  Wilson,  E.L.,  "Numerical  methods  in  finite  element 
analysis,"  Englewood  Cliffs,  N.J.,  Prenti cs-Hal 1 ,  1976. 

43.  Malkus,  D.  S.  and  Hughes,  T.J.R.,  "Mixed  finite  elements  methods  - 
Reduced  and  selective  integration  techniques,"  Computer  Methods  in 
Applied  Mechanical  Engineering,  vol  15,  no.  1,  1978,  pp  63-81 


44.  Lambe  and  Whitman.  Soil  mechanics,  Wiley,  New  York,  1969 


45.  Shen,  et  al ,  Centrifuge  Consolidation  Study  for  Purposes  of 
Plasticity  Theory  Validation,  University  of  California,  Davis,  19f'.5. 


46.  Naval  Facilities  Engineering  Command.  DM  7.  Soil  mechanics, 
foundation  and  earth  structures.  Alexandria,  Va. ,  Mar  1971. 


47.  Reese, 
mont  clay," 
Ill. ,  1969. 


L.C.  and  O'Neill,  M.W. ,  "Field  tests  of  bored  piles 
American  Society  of  Civil  Engineers  Annual  Meeting, 
(Preprint  1008) 


in  Beau- 
Chicago, 


48.  Peterson,  H.  "Application  o  finite  element  method  in  the  analysis 
of  contact  problems,"  in  Proceedings  of  the  International  Conference  on 
Finite  Elements  in  Nonlinear  Solids  and  Structural  Mechanics,  Geilo, 
Norway,  vol  2,  Aug  1977. 

49.  Herrman,  L.R.,  "Nonlinear  finite  element  analysis  of  frictional 
systems,"  in  Proceedings  of  the  International  Conference  on  Finite 
Elements  in  Nonlinear  Solids  and  Structural  Mechanics,  Geilo,  Norway, 
vol  2,  Aug  1977. 


50.  Herrman,  L.R  and  Al¬ 
so  il  system,"  American 
Convention  and  Continuing 
1978.  (Preprint  3125). 


Yassin,  Z.  ,  "Numerical  analysi 
Society  of  Civil  Engineers 
Education  Program,  Pittsburgy, 


s  of  reinforced 
National  Spring 
Pa. ,  Apr  24-28, 


51.  Herrmann,  L.R.,  "Finite  element  analysis  of  contact  problems," 
Engineers,  vol  104,  Journal  of  the  Engineering  Mechanics  Division, 
American  Society  of  Civil  Engineers,  no.  EM5,  Oct  1978. 


52.  Poulos,  H.G.  , 
the  Geotechnical 
Engineers,  vol  107, 


"Cyclic  axial  response 
Engineering  Division, 
no.  GTI ,  Jan  1981. 


of  single  pile,"  Journal  of 
American  Society  of  Civil 


161 


53.  _ .  "Development  of  an  analysis  for  cyclic  axial 

loading  of  piles,"  Third  International  Conference  on  Numerical  Methods 
in  Geomechanics,  Aachen,  Germany,  vol  4,  Rotterdam,  1979. 

54.  Civil  Engineering  Laboratory.  Technical  Note  N-1541:  Coefficient 
of  friction  between  calcareous  sands  and  some  building  materials,  and 
their  significance,  by  P.J.  Valent.  Port  Hueneme,  Calif.,  Jan  1979. 


162 


DISTRIBUTION  LIST 


AF  HO  FREES  Washington  [X 

AF  API 1  DP  1'.  Wnght-Patterson  AFB.  OH 

AFB  AFSCDEEO  (P  Montoul.  Peterson  AFB.  CO,  HO  MAC  DEEM.  Scott  AFB.  IF;  SAMSO  MNND. 
Norton  AFB  (  A 

AFP  SC  DFB.  Tyndall  AFB.  FL;  HO  AFESCTST.  Tyndall  AFB.  FF:  HO  RDC.  Tyndall  AFB.  FI;  HO  1ST. 

Fyndall  AFB.  FF 
AF  HO  FSDOCMS 

ARMY  ARADCOM  SI  INFO  Do.  Dover.  NJ 

ARMS  Comm  Cmd.  lech  Ref  Do.  Huachuca.  AZ;  ERADCOM  Tech  Supp  Dir  (DEFSD-F).  Ft  Monmouth. 

NJ.  FI  SA  EM  ( Kra(evvski).  Ft  Belvoir.  VA:  HODA  (DAEN-ZCM) 

ARMS  CFRI  CF.RF-Z.N.  Champaign.  IF:  library.  Champaign  IF 

ARMS  CORPS  OF  ENGINEERS  HN'DED-C  S.  Huntsville.  AF;  HNDED  SS  .  Huntsville.  AF.  l  ibrary.  Seattle. 
WA 

ARMS  C  RRl  1  l  ibrary.  Hanover.  Nil 

ARMY  EAVES  l  ibrary.  Vicksburg  MS.  WESCV-Z  (Whalin).  Vicksburg.  MS;  WESGP-E  (Green).  Vicksburg. 

MS 

ARMS  ENGR  DISC  Phila.  Fib.  Philadelphia.  PA 

ARMS  MAI  \  MECll  RSCH  CEN  DRXMR-SM  (Fenoel.  Watertown.  MA 
ARMS  MISSIFE  RAD  CMD  Ch.  Docs.  Sci  Info  Ctr.  Arsenal.  At. 

ARMS  BEl.VOIR  RAD  (H  R  STRBE-AALO.  Ft  Belvoir.  VA;  STRBE  BEORE.  Ft  Belvoir.  VA. 

S I  RBI  -C  FIO.  Ft  Belvoir.  VA.  STRBE-WC.  Ft  Belvoir.  VA 
Bl'RFAl  OF  RECLAMATION  D- 1512  (GW  DePuy).  Denver.  CO 
CBC  Feeh  library.  Gulfport.  MS 

(  NO  Code  NOP-OM.  Washington.  DC:  Code  OP-W7J.  Washington.  DC 
COMCBFANI  (  ode  S.VI.  Norfolk.  VA 
COMNAVSl  RF1.ANI  Code  N42A  Norfolk.  VA 
COM(X  EANSYSI.AVI  Fac  Mgntl  Offr.  PWD.  Norfolk.  VA 
COMIRAI  AN  I  SCF.  Norfolk.  VA 

DI  A  DB-hF.  I .  W  ashington.  DC.  DB-6E2.  Washington.  DC:  VP-TPO.  Washington.  DC 

DIRSSP  lech  Fib.  Washington.  DC 

DI  SH  Arms  l  ogistics  Mgt  Center.  Fort  Lee.  VA 

DNA  SITI  IF.  Washington.  DC 

DO!  W  ind  Ocean  Feeh  Div .  Tobacco.  MD 

Dl  1C  Alexandria.  V  A 

DINSRIX  Code  4111.  Bethesda.  MD:  Code  172.  Bethesda.  MD 

FAA  Code  APM-'4il  |  Tortilla  F  Washington.  DC 

I  Ml  I  AN  I  (  EC  Offr.  Norfolk  VA 

FOREST  SERVICE  Engrg  Staff.  Washington.  DC 

GIDEP  OK  Corona.  CA 

GSA  Code  FAIA.  Washington.  DC 

IRE-HTD  Input  Proe  Dir  (R.  Danford).  Eagan.  MN 

LIBRARY  OF  CONGRESS  Sci  A  Feeh  Div.  Washington.  DC 

MARCORPS  FIRST  FSSG.  Fngr  Supp  Offr.  Camp  Pendleton.  CA 

MCDEC  PWO.  Ouantico.  VA 

MCL.B  PWO  (Code  B52D).  Barstow.  C'A 

NAS  Fae  Plan  Br  Mgr  (Code  1S3).  NF  San  Diego.  CA 

NAIF  BCREAC  OF  STANDARDS  R  Chung.  Gaithersburg.  MD 

NAIl  RESEARCH  COUNCIL  Naval  Studies  Board.  Washington.  DC 

N.AVC  OASTSYSCEN  CO.  Panama  City.  FF:  Code  2360.  Panama  City.  FF;  Code  423.  Panama  City.  FF:  Feeh 
Library.  Panama  Cite.  FF 

NAVCONSI'RACEN  (ode  mil' 15.  Port  Hucneme  CA 

NAVFACENGCOM  (ode  03.  Alexandria.  VA:  Code  03T  (Essoglou).  Alexandria.  VA:  Code  04M.  Alexandria. 
VA.  Code  FPO-3C.  Alexandria.  VA:  Code  OOM124  (Feeh  Fib).  Alexandria.  VA:  Code  IIKI2B.  Alexandria. 
VA.  Code  I13C.  Alexandria.  VA 

NAVFACENGCOM  -  CUES  DIV  (  ode  101.  Washington.  DC:  (  ode  405.  Washington.  DC:  Code  406(  . 
Washington.  DC.  Code  407  (D  Schecsele)  Washington.  IX’:  Code  FPO-IC.  Washington.  DC:  Code 
FPO  IPI  .  Washington.  DC 

NA VFACE NG(  ( )M  I  AN  I  DIV  Br  Ofc.  Dir.  Naples.  Italv 

NAIl  BCR t: Al'  OF  STANDARDS  Bldg  Mat  Div  (Rossiterl.  Gaithersburg.  MD 

NAVFACENCiCOM  PANT  DIV  l  ibrary.  Norfolk.  VA 

NAVFACENGCOM  ■  NOR  I'll  DIV  (ode  04.  Philadelphia.  PA  (ode  04AI  .  Philadelphia.  PA.  (ode  40S  Al  . 
Philadelphia.  PA 

NAVF'ACFNGCOM  -  PAC  DIV  (  ode  101  (Kvil.  Pearl  Harbor.  III.  (  ode  0*>P  Pearl  Harbor.  III.  (  ode  2011 
Pe.irl  Harbor.  HI.  (ode  402.  RDIAE  FnO.  Pearl  Harbor.  HI.  labrarv.  Pearl  Harbor.  III 


163 


N  NYI  NCI  NOCOM  -  SOI  III  l)l\  l  ode  hi:,  (  harkslon.  SC  .  (  mk-  405.  Charleston.  SC  .  C  ode  41  »<», 
Ch.iileston.  SC  .  (ieoieeh  Section  (Code  l< C'l.  Charleston.  S<  lihrari.  Ctia rk’ston.  S< ' 

N AN  I  AC  I  NCR  OM  WISI  DIN  C  ode  04K.  Sail  Bruno.  <  \:  I  ihiarv  |(  'ode  04  \:  Jl.  S.m  Kruno.  (  \ 
RDlAl  ln<).  S.m  Kruno.  (  A 
NAN  I  Ul  MICDNI  COMRAl  IS  NRUK'C,  Ouanlieo.  \  \ 

N  woe  I  AM)  Code  cClHi  |M  P.nee).  Ikes  Si  loins.  MS 

\  woe  I  NNSNSCl  N  Code  541  i  ll.tehnian).  S.m  Dieao.  (A.  Code  ll<>4  (leeli  lihi.uvl.  S.m  Dieuo.  (  \.  (ode 
l<(>4:|)  ilknside  I  il'r.m  I.  S.m  Dieeo.  (  A 
S  NVPOSCOl  C  ode  tdWl  (()  Wilson).  Momeiei.  (  A 
\  W  SC  ( )l  C  I  (  ( >1  I  C  ode  C  ">5.  Port  lliieneme.  C  A 

N  MSI  NSNSCOM  Code  0n||4.  N\ .ishmelon.  DC.  Code  CM -ID:.'.  N\ ..sluneton  DC 

N  W  SHIPS  D  Code  :u:  4.  1  one  Ke.ieli.  (  A 

\  \\S|  \  Du,  I  net  Dt\ .  PNNI)  iC  osie  INdooi.  M.uporl.  II 

NNNSlRIWPNClN  (ode  I'll  It  House  l.  D.ihlcren.  \A 

N  \N  NN  NRCOI  I  ue  Coord  iCode  Ml.  Newport.  Rl 

NANNNPNCIN  Code  >50.  Chin.!  lake.  (  A 

\(  )  \  N  Joseph  N  . ulus.  Roekvillc.  Nil) 

NORDN  Head.  Cieoieih  III  (Code  '(''I.  liar  St  lotus.  MS 
NRI  C  ode  5 Si K i .  Washington.  DC 
Nst  (  ode  54  I .  Noitolk.  \  \ 

I’NIIC  Code  IMs.  Point  Mueii.  C  A 

PNNC  C  otle  PM  llilnurO.  Oakland.  ON.  Code  CSC  .  San  Dieito.  C  A.  Code  4:o.  Client  lakes.  II  .  I  iPi.im 
(Code  I  '4 1 .  Peat!  Harbor.  III.  lihrarv.  (iu.mi.  Nlanana  Islands;  lihrarv.  Norlolk.  \  A.  lihrarv,  Pens.ieol.i. 
II  .  lihrarv,  Nokosuka  .IN.  I  eeh  lihrarv.  Sul'ie  II. iv.  RP 
l  S  |)l  PI  Ol  IN  1 1  RIOR  liur  ol  I  and  Menmt  (Code  5S.M.  Washington.  DC;  Nat  I  Park  Sve.  RMR  PC  . 
Denver.  CD 

IN  (dOIOC.K  A I  SI  RNI  N  Marine  (koine)  Oik  ( Pitek  ki )  Reston.  VA 
l  sD  N  I  M  Serv  t  I  Nlaherl.  W  ashmelnn.  DC  .  l  ores]  Sen.  Ree  S.  Atlanta.  O  N 

l  NN  \  (  haitman.  Nleeh  liters;  Dept  Nnnapolis.  Nil).  Met.  I'nere.  Civil  Specs  Hi.  Annapolis.  Ml) 

CM  II  Dl  PI  Ol  NAYIOAMON  A  CHIAN  DIN  O  Amtstrone.  Saeramento.  (A 

CAMIORNIN  SI  Nil  l  N|N  I  RSI  IN  I  netev  leeli  Dept  iKolianl.  Menlo  Park.  CN 

CIIN  ( )(  IINIRMORl  Protect  I  ner  tD.mkins).  I  iverntore,  CN 
Cl  NRKSON  (Oil  Ol  11(11  Cl  Dept  (liaison).  NN 

(OIORNI)O  SCHOOI  Ol  NI|N|  N  Depl  ol  l.nere  (C  Inline ).  (.olden,  CO 

(  Ol  ORADO  SI  N  1 1  CNIYIRNin  (I  Dept  iNelson).  Id  Collins.  CO.  Cl  Depl  iW  C  hat  lie ) .  Port  Collins. 
Nil) 

CORN1I1  I  NIN  I  RSI  I  A  lihrarv.  Ithaca.  NN  .  (  nil  A  I  nviroit  Ineretl  Knllmavi.  Ithaca.  NN 
DNN1IS  A  NIOORI  IIKRNRY  I  os  Nneelcs.  C  A 

1IORIDN  Nil  NNIK  l  N  IN  I  RNl  I  Y  Ocean  I  nere  Depl  I  Nle  Allistet  I  lloea  Raton.  II  .  (ke.m  I  Here  Dept 
tSul.  lloea  Raton.  I  I 

I  I  OKI!)  A  INS  1 1  1 1  1 1  Ol  II  (  IINOI  OOY  (  I  Depllkal.inanl.  Melbourne.  H 

Ol  ORC,|  A  INS  1 1 1 1  II  Ol  lie  IINOI  OON  Null  (  o|  ( Kenton  I.  Ntl.mta  ON  (I  Seol  Naim).  Nll.inla.  ON 
(  I  Seol  i  Ma/anni.  Atlanta,  t  r  N 

[Ns  1 1  I  l  II  Ol  MARINI  S(  II  N(  I  N  lihrarv.  Port  Mamas.  I\ 

IOW  N  SINII  l  N  IN  I  RSI  IN  Cl  Dept  tll.uulvl  Nines.  I  N 
I  NWRINC1  IINIRMORl  N  N  I  I  I  Nil  I  I  lokar/.  I  ivemuuc.  C  N 

I  I  HIOII  IN1NIRSIIN  (I  Dept  llvdranlie'  lab.  Ilethleltem.  PN.  I  inderm.m  1  il't  Set  l  .italoenet . 

Bethlehem.  PA.  Nlarine  Oeoteeh  lab  I  N  Rnhaidsi  Bethlehem.  PN 
NIB  HIC, AN  1 1  (  IINOI  OOIC  Nl  1NINIRSIIN  (I  Dept  Hlaasi  lloiiehton  Nil 

Nil  I  I  net  a  I  ib.  (  ainhtulee  Nl  N.  I  ib.  leel)  Reports.  C  .mihndce.  NI  N  RN  Whitman,  t  .mihudee.  M  N 

NN  CIIN  (  ( 1  Nl  N II  NIIN  COl  |  I  (if  I  il’t.itv  .  Hiooklvn.  NN 

PI  NNSNI  \  NN|  N  SI  Nil  I  N  |N  I  RNI  I  N  OoioW.i  Iniveisi;..  p,il.  p\  Rs.  h  I  I  Snv .!.  i  ■  S'  e.  (  .  II, 

P  N 

PI  RDI  I  ININIRSIIN  (I  Seol  l  Nits.  h.Ktlli  lalavetie  IN  (  I  Seol  1 1  conn.;,  i  Oil. .veils  IN  Ineielih 
I  alavetle.  IN 

SNN  DllC.C)  SI  Nil  ININ  Cl  Dept  i  Ki  ishnaniooi'hv  i  s.m  Di.t"  l  N  (I  I k  pi  <  Noomiiv  t  s.m  Die  c 
SI  N  I  1 1  I  ININIRSIIN  (I  Dept  '  Selnv  ..eelerl  Seattle  WN 

SC)l  i  II W I  SI  RSC  II  INSI  I  neieetis  Svs  Dept  il  -pn/.u  S  ui  Nriioiuo  |\  I  llok.mson  Sm  Nntom.*  IN 
km;  S.m  Nntonio.  I\.  R  DeHart.  S.m  Nntonio  IN 
SINK  ININ  Ol  Nl  W  N’t  )KK  (I  Dept  Hnltalo  NN 

II  \  NS  NANI  l  N|N  I  RSI  I  N  (I  Dept  1 1  edlvttei  I  t  .dfeee  Station.  IN  t  I  I  k  pi  i  Nied/vves  ki  i  t  oik  e. 

Station  IN  (ke.m  I  net  I’roi  (  olke’e  Station  IN 

ININIRSIIN  Ol  (NIIIORNINCI  Dept  Hieivviiki  Belkekv  (Nil  I  )ept  I  Nhu  In  II  i  Ikik.liV  (  \ 

I  rrere  I  W  rib. muon).  Berkekv  t  \  Dans  I  ne'fe  I  )ep!  iDiiik.hu  Bok.I.v  t  N 
ININIRSIIN  Ol  II  I  INOlS  (  I  Depl  ill. dll  lih.nn  II  I  I  I  k  p!'  >  W  ( ..nnhk  i  I  ilvina  II  lib.  nv 

I  rb.uia  II  M  l  Davisson.  I  rh.in.t  II  Nlel/  Rel  llm  I  rb.ma  II 


m 


m 


i 

I 

•S? 

$ 


l  NI\ I  Rsl  I  X  Ol  MASS M  ill  SLITS  Ml  Dept  ( I  leroneumus).  Amherst.  MA 

l  MMRSIH  Ol  MU  HU. AN  I  I  Dept  (Rielurt).  Ann  Arbor.  Ml 

l  NIX  I  Rsl  IX  t)f  N |  HRASKA I  INCOIN  Ross  lee  Shell  l’ro|.  Lincoln.  NI 
I  N|\l  Rs|  I  X  Ol  MAX  Mi  MCI)  NMLRI  (F.ilkl.  Alhtujue  |ue.  NN1 

L  NIX  l  Rsl  I  X  Ol  NOIRi  DAMP  l  I  Dept  lK.nton.il.  Norte  Dame.  IN 

l  N|\ |  Rsl  IX  Ol  RIIODl  ISI  AND  l  !  Dept  (Kovaes).  Kingston.  Rl 

l  N I X  I  Rsl  1  X  Ol  HAAS  XI  XI  SUN  Breen.  Austin.  I  V  Cl!  Dept  (R  Olson).  Austin.  I  \.  I  I  Dept 

,  1  hompsoni.  Xustm.  I  \ 

l  N | X  1  Rsl  I  X  Ol  XX  XSHINlilON  l  l  Dept  iMattoekl.  Seattle.  XX  X.  Rl  lerrel.  Lverelt.  XX  A 

l  N I X  1  Rsl  1  X  Ol  XXISCONSIN  C. reat  lakes  Studies.  <  tr.  Milwaukee.  XXI 

XMI  RU  XN  (  ONI  Rl  II  INS  1 1 1 1  II  Librarv .  Detroit.  Ml 
XMl  II  K  OHS1IORI  RSC  II  Santa  Barhara.  (  A 
XII  XN  l  It  RU  Hill  ID  CO  Rl  Smith.  Dallas.  IX 

XI  SIR  XI  I  X  Swlnex  l  niv  .  Seol  Cl  A  Mine  (I’oulost.  Svdnev 
IIMirill  D  I  rmk.  Columbus.  Oil 

HI  1111  I  HIM  S  1 1  I  I  CO  I  nitre  Dept  iDismukel.  Bethlehem.  PA 
HRIIlslI  I  MB  XSSX  Set  sk  leelt  Dept  (XXilkui'l.  Washington.  DC 
HROXXN  \  ROC  II  Ward.  Houston.  I  X 
l  X N  \ D  A  X  latent  De  Champlain.  D.S  A  .  Matane  Canada 
C  III  X  RON  (HI  III)  RSC  II  CO  Brooks,  la  Habra.  (A 
C  t  >N(  Rl  II  11C  II  C  ORB  A  Anderson,  Laconia.  XX  A 
CONRAD  XSStIC  I  uisom.  X  alt  Nuts.  (A 
C  ( INS  I  Rl  C  I  ION  I  III  I  [  XBS  XL  Liorato.  Skokie.  11 
DRAW)  (ORB  Wright.  Pittsburg.  BA 

I  X  XI  l  XI  ION  XSSOC  INC  M  X  Fedele.  King  of  Prussia.  BA 

II  C.RO  IN  If  R.  INC  l  ibrnrs .  Houston.  IX 

CH  one  UN  IC  XI  [  NCHNi  [  RS  INC  (RF  Murdoekl  Brmeipal.  XVmehestet.  XIX 

( .Rl  XINIXN  XI  R<  ISI’ AC  [  CORP  Teeh  Into  C  tr.  Bethpagc.  N't 

HAI  FA  X  ALDRIC  H.  INC  HP  Aldrich.  Jr.  Cambridge.  XI  X 

Nl  SC  1)1  I  1  tbrars  IC  ode  45X31  Newport.  Rl 

KXXIXN  SC  II  NC  I  S  CORP  D  Saehs.  Santa  Barbara.  CA 

I  XXIONI  DOHI  RIA  CiFOIOCilC  XI  OBSLRX  A  lORX  MeC  oc.  Palisades.  NX 
I  IN  OFISHORI  LNC.RCi  P  Chow.  Sail  Franeiseo  C  A 
1  INI)  X  II  XI  I  I  1HRARX  Doe  Dept  Kansas  Cm.  MO 
X|  XR  X  I  IK  IN  Oil  C  < )  Houston  I  X 

XIOBIl  R  X  D  CORP  Offshore  Lite  librarc.  Dallas.  IX 

XII  I  SI  R  Rl  1111  I  DO|  ,  CONS!  I  TINO  I  NOBS  Riebanls.  New  Xork.  NX 
BROI  s\  (  s  INDISIR1IS  INC  Dir.  Roots  1 1  wins).  Houston.  IX 

Bill  IBs  XSSOC  B  X  Phelps.  Rheem  X'allex.  CA 

BOKII  XNI )  (  I  XII  N|  XSSOC  Cories.  Skokie.  II  .  Rseh  X  Dee  lab  lib  Skokie.  II 

KXXXItiND  IN  1 1  INC  soil  Kelt  Dept  ll  Collet.  Bemtsauken.  NI 

s  XN|  tlX  I  XROKXIORIIS  librars.  I  nermore  C  X 

si  X  II  t  II  t  I  'RB  Ben  mi  Miami  I  I 

SIIXNMIN  X  XXII  sox  INC  libiatiait  Seattle  XXX 

Sill  I  I  Oil  Silt  >RI  INC  I  Doxle  ll.uMoit.  I X 

s|  II  I  I  i'll  c  o  |  x  P  C  iv  il  I  nete.  Houston  IX 

I IDI  XX  X  1 1  R  (  ONS  |  R  c  O  I  I  owlet  X  itemi.i  Be.islt  XX 

I  R  XX  s  X  s  1  |  X I  s  |iu  San  Bel  nai  ditto  C  X  I  net  I  ihr  arv  t  lev  eland  <  >  1 1 

XXI  slINtillol  s|  lllcIRU  CURB  |  ibiatv  Pittsbuie  I'X 

XXIss  I  XNNI  X  IlsINIR  X  XSSOC  I  >  XX  Pkiiei  Notthbiook  11 

XX  l  H  >1 1 XX  XRD  I  I  X  Dl  (  ONsI  I  I  \N|S  R  t  ,..s.  XX  alum  «  reek  (  X  R  I  >■  ■ininene.'  Iloust.u  IN  XX  If 


>.'l 


HRIIXXN  Kill 
HI  I  /I  IX 
Nil  I  U  .If.  i|  i  X 
Nil  K  XII  '  V 


Id  >BI  R  I  I  it i \  v  fsi' 


*.  *.  A  A/.  AV.A  A 


INSTRUCTIONS 


The  Naval  Civil  Engineering  Laboratory  has  revised  its  primary  distribution  lists.  The  bottom  of 
the  mailing  label  has  several  numbers  listed.  These  numbers  correspond  to  numbers  assigned  to  the  list  of 
Subject  Categories.  Numbers  on  the  label  corresponding  to  those  on  the  list  indicate  the  subject  category  and 
type  of  documents  you  are  presently  receiving.  If  you  are  satisfied,  throw  this  card  away  (or  file  it  for  later 
reference). 

If  you  want  to  change  what  you  are  presently  receiving: 

•  Delete  mark  off  number  on  bottom  of  label. 

•  Add  circle  number  on  list 

•  Remove  mv  name  from  all  vour  lists  -  cheek  bos  on  list 

•  Change  my  address  line  out  incorrect  line  and  write  in  correction  (ATTACH  MAILING  LABEL). 

•  Number  ot  copies  should  be  entered  after  the  title  of  the  subject  categories  you  select. 

Fold  on  line  below  and  drop  in  the  mail 

Note  Numbers  on  label  but  not  listed  on  questionnaire  are  for  NCEL  use  only,  please  ignore  them. 
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DISTRIBUTION  QUESTIONNAIRE 

The  Naval  Civil  Engineering  Laboratory  is  revising  its  primary  distribution  lists. 


SUBJECT  CATEGORIES 

t  SHORE  FACILITIES 

2  Construction  methods  and  materials  (including  corrosion 

control,  coatings) 

3  Waterfront  structures  (mamtanance/detanoration  control) 

4  Utilities  (including  power  conditioning) 

5  Explosives  safety 

6  Construction  equipment  and  machinery 

7  Fire  prevention  and  control 

8  Antenna  technology 

9  Structural  analysis  and  design  (including  numerical  and 

computer  techniques) 

TO  Protective  construction  (including  hardened  shelters, 
shock  and  vibration  studies) 

1 1  Soii/rock  mechanics 

13  BEQ 

14  Airfields  and  pavements 

15  ADVANCED  BASE  AND  AMPHIBIOUS  FACILITIES 

16  Base  facilities  (including  shelters,  power  generation,  water  supplies) 

17  Expedient  roads/airfields/bridges 

18  Amphibious  operations  (including  breakwaters,  wave  forces) 

19  Over  the  Beach  operations  (including  containerization, 

materiel  transfer,  lighterage  and  cranes) 

20  POL  storage,  transfer  and  distribution 
24  POLAR  ENGINEERING 

24  Same  as  Advanced  Base  and  Amphibious  Facilities, 
except  limited  to  cold  region  environments 


TYPES  OF  DOCUMENTS 

85  Techdau  Sheets  86  Technical  Reports  and  Technical  Notes 

83  Table  of  Contents  8c  Index  to  TDS 


28  ENERGY/POWER  GENERATION 

29  Thermal  conservation  (thermal  engineering  of  buildings,  HVAC 

systems,  energy  loss  measurement.  pow«r  generation) 

30  Controls  and  alectricai  conservation  (electrical  systems, 

energy  monitoring  and  control  systems) 

31  Fuel  flexibility  (liquid  fuels,  coal  utilization,  energy 

from  solid  waste) 

32  Alternate  energy  source  (geothermal  power,  photovoltaic 

power  systems,  solar  systems,  wind  systems,  energy  storage 
systems) 

33  Site  data  and  systems  integration  (energy  resource  data,  energy 

consumption  data,  integrating  energy  systems) 

34  ENVIRONMENTAL  PROTECTION 

35  Solid  waste  management 

36  Hazardous/toxic  materials  management 

37  Wastewater  management  and  sanitary  engineering 

38  Oil  pollution  removal  and  recovery 

39  Air  pollution 

40  Noise  abatement 

44  OCEAN  ENGINEERING 

45  Seal,  ^or  soils  and  foundations 

46  Seafloor  construction  systems  and  opera  _>ns  (including 

diver  and  manipulator  tools) 

47  Undersea  structures  and  materials 

48  Anchors  and  moorings 

49  Undersea  power  systems,  electromechanical  cables, 

and  connectors 

50  Pressure  vessel  facilities 

51  Physical  environment  (including  site  surveying) 

52  Ocean-based  concrete  structures 

53  Hyperbaric  chambers 

54  Undersea  cable  dynamics 


82  NCEL  Guide  &  Updates 
91  Physical  Security 


□  None- 

rcmovf  my  name 


! 


PLEASE  HELP  US  PUT  THE  ZIP  IN  YOUR 
MAIL!  ADD  YOUR  FOUR  NEW  ZIP  DIGITS 
TO  YOUR  LABEL  (OR  FACSIMILE), 
STAPLE  INSIDE  THIS  SELF-MAILER,  AND 
RETURN  TO  US. 
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